IQQ K I .2708-6437(Online) 
I oi IN *2708-6429(Print) 



Journal of 
Engineering Advancements 


Editor-in-Chief: 

Prof. Dr. Mohammad Mashud 


Volume 01 Issue 03 



SciEn 


Published by: 

SciEn Publishing Group 



Journal of Engineering Advancements 

Apt. #6C-D, House# 191 
Road # 9/A, Dhanmondi R/A 
Dhaka-1209, Bangladesh 

Email: jea@scienpg.com Website: www.scienpg.com/jea/ 


Editor-in-Chief 

Prof. Dr. Mohammad Mashud 

Khulna University of Engineering & Technology 

Khulna-9203, Bangladesh. 

Tel: +880-41-769468 Ext. 405 
Email: mdmashud@me.kuet.ac.bd 


Executive Editor 

Dr. Md. Arifuzzaman 

Khulna University of Engineering & Technology 
Khulna-9203, Bangladesh. 

Tel: +880-41-769468 Ext. 431 
Email: arif48@me.kuet.ac.bd 


Journal of Engineering Advancements Editorial Board Member 


Dr. Abul Mukid Mohammad Mukaddes 
Shahjalal University of Science and Technology 
Email: mukaddesl975@gmail.com 
Bangladesh 

Dr. Chu Chi Ming 
University Malaysia Sabah 
Email: chrischu@ums.edu.my 
Malaysia 

Dr. Mohammad H. Rahman 
University of Wisconsin-Milwaukee 
Email: rahmanmh@uwm.edu 
USA 

Dr. Sivakumar Kumaresan 
University Malaysia Sabah 
Email: shiva@ums.edu.my 
Malaysia 

Dr. Md. Mizanur Rahman 

World University of Bangladesh 

Email: mizanur.rahman @ mte. wub. edu.bd 

Bangladesh 

Dr. Riaz U. Ahmed 

University of Wisconsin-Green Bay 

Email: ahmedm@uwgb.edu 

USA 

Dr. Kazi Mostafijur Rahman 

Khulna University of Engineering & Technology 

Email: mostafij@me.kuet.ac.bd 

Bangladesh 

Dr. Md. Rashedul H. Sarker 
University of Indianapolis 
Email: sarkerm@uindy.edu 
USA 

4 SciEn 


Dr. Seock Sam Kim 
University Malaysia Sabah 
Email: sskim@ums.edu.my 
Malaysia 

Dr. Sabuj Mallik 
University of Derby 
Email: s.mallik@derby.ac.uk 
UK 

Dr. Mohd Suffian Bin Misaran 
University Malaysia Sabah 
Email: suffian@ums.edu.my 
Malaysia 

Dr. Zahir Uddin Ahmed 

Khulna University of Engineering & Technology 

Email: zuahmed@me.kuet.ac.bd 

Bangladesh 

Dr. Mohammad Ilias Inam 

Khulna University of Engineering & Technology 

Email: iliasinam@me.kuet.ac.bd 

Bangladesh 

Dr. Md. Mahfuz Sarwar 
AECOM 

Email: mahfuzsarwar@yahoo.com 
Australia 

Dr. Md. Shariful Islam 

Khulna University of Engineering & Technology 

Email: msislam@me.kuet.ac.bd 

Bangladesh 

Dr. Md. Abdullah A1 Bari 

Khulna University of Engineering & Technology 

Email: abdullahalbari @ me.kuet.ac.bd 

Bangladesh 


Published in: October 2020 
Published by: SciEnPG 


Price: Each Issue BDT 200.00 (US$ 15) 


ISSN: 2708-6437 (Online) 2708-6429 (Print) 




Journal of 

Engineering Advancements 


Editor-in-Chief 

Prof. Dr. Mohammad Mashud 

Department of Mechanical Engineering, 

Khulna University of Engineering & Technology, Khulna, Bangladesh 


Executive Editor 

Dr. Md. Arifuzzaman 

Department of Mechanical Engineering, 

Khulna University of Engineering & Technology, Khulna, Bangladesh 



SciEn 


Published by: SciEn Publishing Group 

Apt. # 6 C-D, House # 191, Road # 9/A 
Dhanmondi, Dhaka-1209, Bangladesh 
Email Address: jea@scienpg.com 

www.scienpg. com/j ea/ 


This page is intentionally left blank 



Journal of Engineering Advancements 


Volume 01, Issue 03 October 2020 

CONTENTS 

Original Articles 

01. Post-consumer PET Bottle Recycling: Chemical Dose Optimization 

Md. Shahadat Hossain and Md Anisur Rahman . 70 

02. Finding the Reasons for the Delay Time in a Highway by Analyzing the Travel 
Time, Delay Time and Traffic Flow Data 

T M Junaid Bashar, Md. Sabbir Hossain and Shah Istiaque . 76 

03. A Numerical Approach to Investigate the Influence of Resonator Setting and 
Volume Fraction on Stop Bands in an Acoustic Metamaterial 

Riaz Ahmed, Hossain Ahmed and Sourav Banerjee . 85 

04. Sustainability Analysis of Different Types of Power Plants Using Multi-Criteria 
Decision Analysis Methods 

Anik Saha Dipto, Md Abdullah Al Bari and Shadman Tahsin Nabil . 94 

05. Combustion and Emission Characteristics of a Diesel Engine Operating with 
Varying Equivalence Ratio and Compression Ratio - A CFD Simulation 

Kazi Mostafijur Rahman and Zobair Ahmed . 101 


iii 








This page is intentionally left blank. 



https://doi.org/10.38032/jea.2Q20.Q3.QQl 


4 SciEn 

Journal of Engineering Advancements Vol. 01(03) 2020, pp 70-75 


Post-consumer PET Bottle Recycling: Chemical Dose Optimization 

Md. Shahadat Hossain* and Md Anisur Rahman 

Department of Chemical Engineering and Polymer Science, Shahjalal University of Science and Technology, Sylhet - 3114, 

Bangladesh 

Received: August 17, 2020, Revised: September 11, 2020, Accepted: September 12, 2020, Available Online: September 15, 2020 


ABSTRACT 


Polyethylene terephthalate (PET) bottles are being used in our daily life and consequently go to the landfills after their use. Additionally, 
virgin PET resins are produced from nonrenewable resources, such as fossil fuels, whose reserves are depleting continuously. Therefore, 
to maintain ecological and environmental balances as well as for sustainable development, post-consumer PET (pcrPET) bottles should 
be recycled. Among many recycling processes, mechanical recycling of pcrPET is attractive due to lower cost involvement. One of the 
most crucial and important processes of mechanical recycling is hot washing for contaminants removal. Hot washing uses a cleaning 
solution made of caustic soda (NaOH) and detergent at elevated temperature. In this paper, caustic soda and FORYL LFO (FLO) 
detergent doses were changed gradually to investigate effective contaminants removal through colorimetric study. Finally, 
concentration vs. absorbance graphs from colorimetric study suggests that 2% NaOH and 2% FLO detergent is the optimum chemical 
dose at hot washing for pcrPET recycling. 


Keywords: PET Bottle; Recycling; Chemical Dose; NaOH; Detergent; Colorimetric Analysis. 
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1. Introduction 

Plastic bottles are widely used at beverage, water, oil, and 
pharmaceutical industries in Bangladesh [1], [2]. Most of them 
are directly goes to the landfill after their first hand use but those 
bottles are not biodegradable [3]. As a result, they create 
potential hazard for our agricultural lands, sewage systems as 
well as for water networks across the country [4]. On the 
contrary, Bangladesh is currently producing different plastic 
products from both imported virgin polymer and locally recycled 
plastics [3], [5]. Bangladesh used 75,000 metric tons (MT) 
plastic raw materials mainly for packaging materials, household 
accessories and producing construction materials in 2015 [5], [6]. 
At the same time, Bangladesh recycled around 65 percent of its 
post consumer plastic materials and this trend was increasing in 
every year [2], [6]. This post consumer plastic recycling provides 
Bangladesh two fold advantages: firstly for plastic waste 
management and secondly for saving money which otherwise 
would have been used for importing virgin polymer [4], [6],[7]. 
Unfortunately, up to now most of the recycled plastic materials 
are used only by small and medium-sized enterprises (SMEs) and 
yet no big plastic conglomerates have used recycled plastic 
materials due to its inferior quality [5]. In a book, United Nations 
ESCAP [6] showed that SMEs used 60 percent of the total 
recycled plastic materials and subsequently produced low quality 
plastic products. Rest of the unused and recycled plastic 
materials was exported to China. For instance, Bangladesh 
exported USD 6 million worth recycled plastic materials 
between July 2008 and January 2010 due to lack of proper 
recycling technology [5], [6]. However, Bangladesh has huge 
demand for quality recycled plastic materials inside the country 
[2]. Therefore, recycled plastic materials quality should be 
improved through proper plastic collection method and proper 
washing for decontamination. 


Among various plastic materials, most widely used plastics 
in Bangladesh are polyethylene (PE), polyvinylchloride (PVC), 
polypropylene (PP), polyethylene terephthalate (PET) and 
polystyrene (PS) [2], [3]. Among which post consumer PET 
(pcrPET) is recycled mostly and has been used by SMEs for 
bottle production for beverage, oil, water, and medicine industry 
as well as for PET sheet and polyester fiber production. But still 
now pcrPET is recycled in Bangladesh by conventional method 
where it is manually sorted firstly, later grounded and finally cold 
washed by detergent and subsequently dried for 
decontamination. But, this conventional recycling method 
produces inferior quality recycled PET flakes in Bangladesh. 

However, it is widely known that during pcrPET materials 
recycling, caustic soda (NaOH) and detergent are used in hot 
washing, but dosage of chemicals is different in different 
countries [8]-[10]. This is because pcrPET bottles are collected 
in wide variety of methods in different counties, which make 
them contain different amounts of dirt and contaminants [lll- 
ll 3]. On the basis of contaminants adherence to pcrPET bottles 
after their first hand use, 1-3 % NaOH along with small amount 
of detergent are used for contaminants removal [8], [10]. Some 
authors also report usage of antifoaming agent in hot washing for 
foaming phenomena reduction [8]. Therefore, if pcrPET bottles 
are recycled in Bangladesh, chemical dosage in hot washing 
should be determined by considering amount and types of 
contaminants present in those post consumer bottles [14], [15]. 
If improper amount of washing chemicals are used in hot 
washing step of recycling, it will give yellowish color in 
subsequent recycling steps, for example, extrusion, injection 
molding, and blow molding due to presence of trace amount of 
NaOH in recycled PET materials [13], [16]-[19]. Therefore, 
determination of chemical dosage in hot washing is crucial for 
total pcrPET bottle recycling. Use of proper chemical dosage will 
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turn potential environmental hazard into valuable PET flakes and 
subsequently will reduce the environmental pollution through 
PET bottles recycling. In this paper, pcrPET bottles generated 
across the country were collected and then a step by step 
recycling process was developed for contaminants removal by 
employing chemical dosage. Different chemical dose - varying 
proportion of NaOH and FORYL LFO (FLO) detergent - were 
tested in hot washing, a crucial step of recycling process, to find 
out the optimum chemical dose for pcrPET bottles recycling. 

2. Methodology 

2.1 Experimental Set-up 

It was assumed in this study that a pcrPET recycling plant 
is running with a production capacity of 1.0 MT recycled PET 
flakes per hour where 30 minutes of hot washing for pcrPET 
flakes is employed. A hot washer (is shown in Fig. 1) with a 
100 rpm stirrer uses 6 m 3 water and NaOH plus FLO detergent 
dose for pcrPET flakes hot washing. Therefore, it is estimated 
that every 30 minutes hot washer contains 0.50 MT flakes 
within 6 m 3 of water for hot washing. However, this paper used 
small scale hot washer accommodating 1 liter (L) water for 30 
minutes of hot washing. Consequently, 0.083 kg pcrPET flakes 
were used for hot washing with NaOH and FLO detergent 
dosage at elevated temperature. 



Fig. 1 Schematic of a hot washer. 


2.2 Sample Collection 

pcrPET bottles sample was collected from five big suppliers 
who collected bottles through their small collection hub across 
the country and 2 kg sample was collected from each supplier. 
Moreover, another 2 kg pcrPET bottles sample was collected 
directly from waste collection bins and landfills. Later, all the 
samples, i.e., 12 kg bottles were mixed together to homogenize 
pcrPET bottles generated throughout the country. However, 
each time 0.1 kg pcrPET bottles were taken from homogenized 
sample for subsequent sample pretreatment. 

2.3 Sample Pretreatment 

pcrPET bottles taken from homogenized sample also contained 
other plastics such polyethylene (PE), polypropylene (PP), 
polyvinyl chloride (PVC) etc. and metals alongside pcrPET due 


to lack of proper post consumer plastic waste collection method 
in Bangladesh. The undesired plastics were removed manually 
in this step. 

2.4 Cold Washing 

After initial pretreatment, pcrPET sample were washed 
with cold water at ambient condition to remove dirt from the 
sample in a rector called cold washer. Cold washer was mainly 
a 2 L bucket containing a slow speed stirrer (50 rpm) which 
accommodated initially taken pcrPET sample. Cold washing 
provided relatively clean pcrPET bottles which later chipped 
into desired size in subsequent processing steps. 

2.5 Sample Cutting 

Cleaned PET bottles sample were chipped into 7-13 mm 
sized pcrPET flakes (Fig. 2) in accordance with the recycled 
PET flakes size of several large recycled PET flakes 
manufacturer from Asia to Europe [20]-[22]. Sample cutting 
was carried out without any kind of drying. 



Fig. 2 pcrPET flakes sample for hot washing after cutting into 
desired 7-13 mm size. 

2.6 Gravity Separation 

Gravity separation in another water reactor was employed 
to separate PE, PP, PVC etc. materials from pcrPET flakes 
which finally provided flakes sample for hot washing. This 
flakes sample contained mostly glue and dirt contaminants 
which were removed in hot washing. 

2.7 Hot Washing 

Hot washing used NaOH and FLO detergent dose at 80 °C to 
remove contaminants from pcrPET sample flakes while lower 
and higher temperature was avoided due to ineffective washing 
and excessive foaming correspondingly. Hot washing was 
carried out in small scale reactor which has 30 minutes 
residence time for sample flakes. Moreover, this small reactor 
also contained a high speed stirrer (100 rpm) for effective 
contact between flakes and hot washing chemicals. Chemical 
dose - NaOH and FLO detergent percentage - was changed 
gradually in both reactors until hot washing produced 
contaminants free quality recycled PET flakes. 
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Actually, NaOH and FLO detergent dose percentages were 
employed as hot washer reactive volume (1 L). For example, 
detergent dose was slowly augmented from 1% to 3% percent 
of total hot washer volume (1 L) but NaOH dose was constant. 
Three different constant NaOH dose values such as 1, 2 and 3% 
were selected. Each time different detergent dose, varying from 
1% to 3%, at a constant NaOH dose was used in hot washing. 
Then optimum chemical dose in hot washing reactor at 80 °C 
temperature was determined from recycled PET flakes 
colorimetric analysis. 

2.8 Sample PET Flakes Washing and Drying 

After hot washing and before characterization by 
colorimetric analysis, recycled PET flakes were washed 
vigorously and dried. This is because the presence of trace 
amount of unused hot washing chemicals (NaOH and FLO 
detergent) can produce yellowish color in recycled PET flakes. 
Therefore, extensive cold washing, for example, centrifugal 
washing, float washing and spray washing were carried out one 
after another for at least 10 minutes to remove unused hot 
washing chemicals. After hot washing chemicals removal, PET 
flakes were dried for 10 minutes at 130 °C for water removal. 
Drying of PET flakes produced contaminants free, dry and 
quality recycled PET flakes which were again dried at 220 °C 
for another half an hour. 

2.9 Recycled PET Flakes Colorimetric Analysis 

After consecutive dryings, recycled PET flakes from every 
chemical dose were stored individually. Stored recycled flakes 
from every single dose were used to prepare pcrPET solution 
with 0.025, 0.05, 0.1, 0.15, and 0.2 concentrations. Then a 
colorimeter was used for absorbance measurement against 
every concentration. Concentration vs. absorbance curve for 
different chemical doses were compared with pcrPET (without 
using any chemical dose) and pure PET to find out optimum 
chemical dose at hot washing. Before that two calibration 
curves, one for pure PET and another one for pcrPET with any 
chemical dose, were prepared for aforementioned comparative 
study. 

3. Results and Discussion 

3.1 Calibration Curve of Pure PET and pcrPET 

Colorimeter absorbance was measured against different 
pure PET concentrations and results are represented in Fig. 3. 
Pure PET means those PET that have not used in any packaging 
or other applications. It is observed that absorbance is gradually 
increasing with the increment of PET concentrations. Results 
of absorbance at different concentrations of pcrPET are also 
shown in the same figure. Similar trend was observed in this 
case as like as that obtained in case of pure PET. However, 
higher absorbance value is achieved in case of pcrPET than 
pure PET as expected. Higher absorbance value is resulted 
from higher contaminations adherence after pcrPET bottles 
first hand use. Furthermore, every concentration of both types 
of PET is also compared. It was observed that for every single 
concentration, pure PET’s absorbance is lower than that of 
pcrPET. This scenario was obtained since pure PET had no 
chance to adhere any kind of contaminations such as dirt, glue, 
PVC etc. as compared to pcrPET. 


1 Pure PET □pcrPET 



Concentration (mol/L) 


Fig. 3 Comparison of colorimeter absorbance between pure 
and pcrPET. 

3.2 Hot Washing with Varying Chemical Dose 

3.2.1 Hot washing with 1% NaOH and 1%, 2%, and 3% FLO 
(detergent) 

Initially, PET flakes hot washing started with lower NaOH 
dose, 1%, but detergent dose was increased gradually from 1% 
to 3%. Among them lowest chemical dose, 1% NaOH and 1% 
detergent, produced blackish recycled PET flakes (Fig. 4 a) 
after consequent drying at 130 °C and 220 °C. NaOH used in 
the hot washing acted as aggressive cleaning agent for pcrPET 
flakes surface but such a lower dose - 1% NaOH - could not 
remove mainly the dirt and glue contaminants from the flakes 
[8], [23]. As a result, discoloration of recycled PET flakes 
happened after consecutive drying since glue and dirt were 
identified as contaminants responsible for flakes discoloration 
by several authors [8], [24] However, blackish color from dirt 
contaminants offsets the slight yellowish color resulted from 
glue contaminants; resulted blackish recycled PET flakes at 
lower chemical dose [25]. Overall, as the lowest chemical dose 
did not remove contaminants such as dirt, glue etc. properly, 
colorimeter analysis of recycled PET flakes was not performed 
at this chemical dose. 

Later, hot washing chemical dose was increased further to 
1% NaOH and 2% FLO level. Recycled PET flakes produced 
through this chemical dose are presented in Fig. 4 b whereas its 
colorimetric analysis is shown in Fig. 5. In Fig. 5, a 
comparative analysis between this recycled flakes with pure 
and pcrPET is drawn. Graphical representation at Fig. 5 shows 
that recycled flakes with this chemical dose stays in between 
the pure and pcrPET flakes. Actually, 1% NaOH and 2% 
detergent dose (1% vs 2%) lowered the absorbance of recycled 
flakes highest 57.89% (at 0.01 mol/L) compared to pcrPET and 
on average this chemical dose lowered 40.25% absorbance. 
This implies that above chemical dose can significantly remove 
contaminations from pcrPET. 

On the contrary, still there is a chance to improve 
contaminants removal as recycled PET flakes line stays above 
the pure PET line as shown in Fig. 5. Therefore, chemical dose 
again increased at 1% NaOH and 3% FLO level. Although 
recycled flakes of this chemical dose were cleared 
satisfactorily, heavy foaming problem was encountered during 
PET flakes hot washing with this chemical dose. Therefore, 
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colorimetric analysis was skipped for this chemical dose too 
since commercial application of this chemical dose is 
impractical due to heavy foaming at hot washing reactor. 



(a) 


3.2.2 Hot washing with 2% NaOH and 1%, 2%, and 3% FLO 
(detergent) 

As Fig. 5 suggests that further improvement is possible for 
contaminants removal, chemical dose was changed, especially, 
NaOH dose was increased to 2% and FLO dose was remained 
same as previous 1%, 2% and 3%. However, lowest chemical 
dose in this combination again produced blackish recycled PET 
flakes and therefore colorimetric analysis was excluded for this 
dose. 

Then subsequent chemical dose combination, 2% NaOH 
and 2% FLO, produced comparatively cleaner recycled PET 
flakes and their colorimetric representation is shown in Fig. 6 
(2% vs. 2%). Fig. 6 shows that highest 62.5% absorbance 
reduction (at 0.05 mol/L concentration) is possible for recycled 
PET flakes with 2% NaOH and 2% FLO dose and on average 
50.29% absorbance reduction is possible compared to pcrPET. 
This means that contaminants removal has been improved 
further due to chemical dose increase. As a result, recycled 
flakes line approaches further closer to pure PET line and away 
from pcrPET line in Fig. 6. 



(b) 

Fig. 4 Recycled PET flakes after hot washing and drying with 
(a) 1% NaOH and 1% detergent (b) 1% NaOH and 2% 
detergent. 



Fig. 5 Colorimetric analysis of recycled PET flakes with 1 % 
NaOH and 2% percent FLO dose and comparative study 
between recycled and pure and pcrPET. 


However, 2% NaOH and 3% FLO produces slightly 
yellowish recycled PET flakes instead of white and cleaned 
flakes after drying at 130 °C and 220 °C. Moreover, huge 
foaming problem was also faced with this chemical dosage. 
Naturally detergents are prone to foam formation but in this 
chemical dose heavy foaming was formed due to addition of 
higher proportion of detergent to the very close to bubbling 
point of water during hot washing. Simultaneously large 
percentage of detergent and proximity of hot washing 
temperature (80 °C) to bubbling point of water were 
responsible for heavy foaming [18], [26]. Due to above 
mentioned problems, commercial application of this chemical 
dose for pcrPET recycling is not suggested and the colorimetric 
analysis of this recycled PET is thus excluded. 



Fig. 6 Colorimetric analysis of recycled PET flakes with 2% 
NaOH and 2% percent FLO dose and comparative study 
between recycled and pure and pcrPET. 

3.2.3 Hot washing with 3% NaOH and 1%, 2%, and 3% FLO 
(detergent) 

Previous hot washings were carried out at a constant 
NaOH dose either at 1% or 2% and varying percentages of FLO 
dose such 1%, 2% and 3%. From previous result it was found 
that blackish recycled flakes were produced at lower chemical 
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dose (1% FLO) and foaming problem was encountered at 
higher chemical dose (3% FLO). As a result, above two 
chemical doses were excluded while NaOH dose was increased 
further to 3%. 

However, 3% NaOH and 2% percent FLO detergent 
produces recycled PET flakes which turned completely 
yellowish after subsequent drying at previously mentioned 
temperatures. This indicates that the extra addition of chemical 
during hot washing was not removed in several washing steps. 
As a result, remaining NaOH on PET flakes surface turned 
yellowish after successive drying at 130 °C and 220 °C. This 
yellowish phenomenon can also be explained by natural 
tendency of NaOH; turns yellowish after drying at elevated 
temperature [10], [19]. In addition, higher usage of NaOH 
incurs higher operating cost to PET recycling. As this chemical 
dose produces undesired recycled PET flakes, colorimetric 
analysis of this dosage has not also been included in this paper. 

3.3 Hot washing chemical dose optimization 

Recycled PET flakes produced through two best chemical 
doses (1% NaOH and 2% FLO, 2% NaOH and 2% FLO) at hot 
washing are compared together with pure and pcrPET in Fig. 7 
to find out optimum chemical dose for hot washing. It can be 
observed that the absorbance of recycled flakes were always 
lower at every concentration compared to that of pcrPET. This 
means that the above two chemical doses are capable of 
removing attached contaminants from pcrPET. However, 2% 
NaOH and 2% FLO dose has lower absorbance value than other 
chemical dose which implies that this chemical dose has higher 
efficiency for contaminants removal from pcrPET at every 
concentration. Moreover, the absorbance of recycled PET 
flakes corresponding to this chemical dose has closest 
proximity to pure PET, which clearly indicates that recycled 
PET flakes produced through this dose is closest to pure PET- 
in terms of quality-than other recycled PET flakes. Therefore, 
2% NaOH and 2% FLO detergent dose is the best suited 
candidate for pcrPET flakes hot washing. 


□Pure PET □ pcrPET ni%vs.2% □2%vs. 2% 



0.05 0.10 

Concentration (mol/L) 


Fig. 7 Comparative study between pure and pcrPET with 
recycled PET at two chemical doses. 

4. Conclusion 

Mechanical recycling is an attractive process across the 
world for pcrPET bottles recycling due to its lower capital and 
operating cost involvement. Successful recycling mainly 
depends of contaminants removal from post consumer bottles 


after their first time use. This study used hot washing for 
effective contaminants removal. For hot washing, NaOH and 
detergent FLO was used in different proportions. Lower 
detergent dose, for example, 1% FLO, could not remove the 
dirt contaminants properly from the pcrPET resulted blackish 
recycled PET flakes. On the contrary, higher detergent dose - 
3% FLO - created heavy foaming during the hot washing 
which made difficult to finish the hot washing and subsequent 
collection of recycled PET flakes for further colorimetric 
analysis. Similarly, higher NaOH dose such as 3% NaOH 
produced yellowish recycled PET flakes (undesired quality) 
after hot washing and further drying. It occurred due to the 
presence of higher amount of NaOH at PET flakes surface after 
hot washing and later extensive cold washing. Finally, 1% 
NaOH and 2% detergent dose and 2% NaOH and 2% detergent 
dose produced desired quality recycled PET flakes. Among 
them 2% NaOH and 2% FLO detergent (a non ionic detergent) 
was found to be optimum chemical dose for hot washing since 
this chemical dose removed the highest percentage of 
contaminants from pcrPET. However, this paper studies only 
chemical dose optimization at the hot washing, but there are 
several other factors that should be considered for complete 
mechanical recycling of post consumer PET bottles. 
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ABSTRACT 


The objectives of this study are to show a comparison among travel time, running time, delay in peak and off-peak hours on different 
days of a week, and reasons behind the delay time. Moving car observer method has been carried out to count the traffic flow, journey 
time, running time, and delay time. Total vehicle flow, and comparative vehicle flow during the peak hour and off-peak hour for 
workdays and weekend days were surveyed to show a relationship between delay time and traffic flow. As the traffic flow increases 
the delay time also will increase. To measure the reasons behind the delay time of Fulbarigate-Daulatpur road, spot speed study was 
done in two intersections of the road. The 15 th percentile speed for Religate intersection road is 10 K.P.H. That means 85% of vehicles 
tends to go faster than 10 K.P.H. in this section. And hence vehicles with less than or equal 10 K.P.H. are responsible for increasing 
the delay time in this area. This speed limit should be prohibited in this section to reduce delay time and congestion. For the Mohosin 
More road intersection, the prohibited speed is also 10 K.P.H. Easy bike and Mahindra account for the congestion of these intersections 
and tend to stop in these intersections to collect passengers which creates unwanted queue in this study area. Controlling traffic flow at 
intersections can be a possible way to reduce the congestion rate of Fulbarigate-Daulatpur road. 


Keywords: Spot Speed; Journey Time; Travel Time; Percentile; Congestion; Delay Time. 
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1. Introduction 

The provision of adequate and equitable service for all 
groups is the main target of urban planning. Transportation is one 
of the most important elements of such a service [1]. 
Transportation is a process that involves the movement of 
commuters, goods, and services from a given point of origin to a 
specific destination [2]. It determines the regional patterns of 
development, economic viability, environmental impacts, and 
maintenance of socially acceptable levels of quality of life. It is 
a means to access business activities, education, employment, 
and recreational opportunities [3]. Various transportation modes 
are available like road-based, railway based, air-based, and 
water-based [4]. In our country, road-based transportation is 
dominant among them. In a road-based transportation system, 
the public transport system is most popular [5]. 

“Since the end of the Second World War, high demand for 
mobility increases due to rapid urbanization and economic 
growth in third world cities” [6]. Like other developing countries 
in Bangladesh, problems related to the lack of traffic mobility are 
getting extremely high [7]. Bangladesh cannot bear the cost of 
the financial and ecological misfortune that came about because 
of this serious traffic obstacle [8]. Main urban areas like Dhaka 
(eats up to 3.2 million working hours/day), Chittagong, and 
Khulna are under the disturbing worry of this issue [9]. Khulna, 
the third biggest city of Bangladesh, in such a case is close to 
repeating a similar incidence [10]. 

In respect of Khulna division’s transportation system, 
Khulna-Jessore Highway is the most important road network for 
both intra and intercity connectivity. “Besides after the 
completion of Padma Bridge transportation between the Khulna 
region and the Dhaka region will become even easier which will 
increase the economic activities of Khulna. To support this 


economic growth and future public demand a well-planned 
transportation system is a prerequisite” [11]. The common modes 
of public transportation in Khulna are pedal rickshaw, baby taxi, 
battery bike and auto-rickshaws etc. [12]. 

However, a comparison among travel time, running time, 
delay in peak, and off-peak hours on different days of a week is 
shown in this study, and the reasons behind the delay time are 
found out. Several traffic studies e.g. volume study, speed study, 
etc. have been carried out for this study. By conducting travel 
time, speed, and delay studies congestion can be measured [13]. 
Moving car observer method is used in measuring travel time, 
delay time, vehicle flow. Delay can be measured by doing field 
measurement, simulation, analytical derivation, or by using a 
combination of these methods [14]. The reasons behind the delay 
time at two intersections of Fulbarigate-Daulatpur road is carried 
out by using spot speed. This type of evaluation of the public 
transport system is very suitable for a highly populated country 
like Bangladesh [15]. As the western part of Khulna is facing 
acute traffic congestion due to the haphazard movement of 
traffic, the amount of delay time has increased. Particularly the 
delay time at major intersections and road sections in Khulna 
Metropolitan City is increasing faster [16]. A huge number of 
slow-moving vehicles create a huge delay time due to the lack of 
proper management [7]. 

There are 66 National highway routes in Bangladesh and 9 
of them are from Khulna division. The route chosen for this 
research is called Khulna-Jessore road and our study area 
Fulbarigate to Daulatpur is within Khulna City Corporation. The 
road covers 3.2 kilometers distance [17]. Thus, this research will 
assist to develop a proper transportation system. In the end, the 
output of this research can be useful in further tasks associated 
with the Khulna-Jessore highway. 
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2. Methodology 

To conduct the study, traffic flow, travel time, running time, 
and spot speed are studied. The study has been conducted based 
on primary data. Primary data sources are - Reconnaissance 
survey, Spot Speed Survey, Running Speed, Journey Speed 
Study, and Delay Studies. Moving car observer method has been 
used for the measurements of running speed and journey speed. 
Congestion and delay time in two intersections of Fulbarigate- 
Daulatpur road is carried out by using spot speed. Spot speed has 
been surveyed through a direct-timing procedure. 

2.1 The moving car observer method 

The field survey was carried out to collect the data on 
Fulbarigate to Daulatpur which has a 3.9 km length. To count the 
traffic flow, average journey Speed, average running speed, and 
delay time, moving car observer method has been used. Several 
test runs are made along the road section and a group of observers 
records various parameters. Four observers were fixed for 
recording four types of data. The four types of data were - 

• Total journey time, delay points and delay time 

• Number of overtaking vehicles 

• Number of overtaken vehicles 

• Number of modes coming from the opposite direction 

The average journey time (min.) for all vehicles in a traffic 
stream in the direction of flow (q) is given by following equation. 

q = (Ma+M w )/(T a +T w ) 


Where, T avg = Average journey time, Tw = Average journey 
time when vehicle is travelling in the stream, 

M w = Overtaking vehicle minus Overtaken vehicle, M a = 
Average number of vehicle in opposite direction, 

Ta = Average journey time when vehicle is travelling in 
opposite stream [18]. 

2.2 Spot speed study 

For our study area, the stopwatch method is preferred. For 
calculating spot speed 60/100 m length at each link of a highway 
will be taken which has free flow. When the vehicles passed at 
the first reference point, the first observer signaled the second 
observer and then the second observer started the stopwatch, and 
finally when the vehicles crossed the other reference points, then 
stop the watch and recorded the time. The speed of a particular 
vehicle will be calculated from how much time it takes to cross 
a certain distance. Thus, from known distance and time, found 
out the vehicle’s spot speed. The spot speed of the vehicles was 
measured four times in two days. 

The cumulative speed distribution curve was prepared to 
analyze the spot speed condition for different speed percentile 
groups. The 15 th , 50 th , 85 th ’ and 98 th percentile speeds are four 
parameters that are commonly used in traffic safety and traffic 
engineering. For example, the 15 th percentile speed reflects the 
minimum or critical speed limits, 50 th percentile speed reflects 
median speed limit, 85 th percentile speed reflects the optimum 
speed limit and 98 th percentile speed reflects the maximum or 
design speed limit [19], [20]. 


Tavg — Ta~(Mw/(j) 


Table 1 Moving car observer data collection sample. 


Morning Weekend Peak Hour 

Focation 

Direction 

Journey Time 
(minutes) 

Stop Time 
(minutes) 

Number of Vehicles 

Over Taking 

Over Taken 

In Opposite 
Direction 

Fulbarigate to 
Daulatpur 

N-S 

6.7 

0.6 

51 

47 

159 

5.8 

0.3 

63 

37 

127 

6.8 

0.5 

55 

51 

161 

6 

0.4 

57 

42 

147 

5.6 

0.2 

68 

39 

138 

6.9 

0.35 

52 

49 

162 

Total 

37.8 

2.35 

336 

265 

894 

Mean 

6.3 

0.4 

56 

44.2 

149 


Source: Field Survey 


3. Data collection & calculation 

As shown in Table 1 data collection sample and sample 
calculation of moving car observer Survey. Same procedure 
carried out both workdays and weekends during Peak Hours 
(8.00-9.00 AM, 12.00-1.00 PM, 5.00-6.00 PM) and Off-Peak 
Hours (10.00-11.00 AM, 3.00-4.00 PM) in both directions. 

4. Data Analysis 

For finding out the possible causes of delay time and the 
level of service of the Fulbari to Daulatpur highway- this study 
focused on several factors. First of all, through the vehicle flow 
survey the amount of total vehicle flow during the peak hour and 


off-peak hour for workdays and weekend days were determined 
to find out if there was any anomaly in total vehicle flow number 
in this road. Secondly using the spot speed survey in different 
intersections in the road the minimum and maximum speed for 
different types of vehicles available was found out. After 
analyzing the data the ideal speed for every vehicle has been 
determined. It also shows those vehicles which travel below the 
15 th percentile speed that causes unwanted congestion in that 
road. Here, to visualize the data of spot- the related speed map 
for the two different intersections is also prepared. Through this 
analysis, the connection between vehicle flow and vehicle speed 
is figured out which eventually leads to those factors which 
impact the delay time in this road. 
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Table 2 Sample calculation of moving observer survey. 


Calculation 

M w = Overtaking vehicle minus overtaken 

11.8 

vehicle 

M a = Average number of vehicle in 

149 

opposite direction 

Tw = Avg. journey time when vehicle is 

6.3 min 

travelling in the stream 

Ta = Avg. journey time when vehicle is 

6.4 min 

travelling in opposite stream 

(Flow) Q veh./hr = 

759 veh./hr 

Avg. Journey time 

5.37 min 

Avg. journey speed = ( 3.9*60) / 5.7 

43.57 Kmph 

Avg. Running time = Avg. journey time - 

4.98 min. 

Avg. Stopped Delay 

Avg. Running Speed = (3.9 * 60)/4.98 

47 Kmph 


4.1 Journey time and delay time study 

This study considers more than one intersection to gather 
data and both of the intersections don’t have any traffic signal. 
That means the delay time found in this study is caused by other 
factors than by the traffic signals. Here the average journey time, 
running time, and delay time is calculated for workday and 
weekend during the pick and off-pick hours. This study considers 
both sides of the vehicle flow (North to South and South to 
North) to analyze the travel time and delay time conditions for 
these intersections. 

As Journey time = Running time + Delay time, the journey 
time is always significantly more than the running time. The 
more the difference between the journey time and running time 
the more the delay time is. Now delay time can be caused by 
many factors including traffic signal, road condition (Ex: 
Potholes, speed breaker, road width), vehicle type (Ex: Car, Bus, 
Auto rickshaw, Rickshaw), total number of vehicles, etc. Now 
this study area doesn’t have any traffic signals that means, other 
factors that cause traffic jams or traffic delays are the reason 
behind the delay time. 


Travel time & delay of PeakHour (N-S) 
Trip 



j- Workday Weekend 


■ Ave rage Journey Time 1 Ave rage Running Time 

■ Ave rage Delay time 

Fig. 1 Travel time & delay of peak hour (N-S) trip. 

As shown in Fig. 1 travel time & delay of peak hour (N-S) 
trip workday has the highest travel time & delay between 
workday and weekend, same for the return trip (S-N) Peak Hour 
workday has the highest travel time & delay as shown in Fig. 2, 


for Travel time & delay of off-peak period (N-S) trip workday 
has the highest Travel time & delay is shown in Fig. 3, Fig. 4 
shows the return of off-peak period trip (S-N). 

Travel time & delay of PeakHour (S-N) 

Trip 



■ Ave rage Journey Time ■ Ave rage Running Time 

■ Average Del ay time 


Fig. 2 Travel time & delay of peak hour (S-N) trip. 

Travel time & delay of Off-Peak Hour 
(N-S) Trip 

« 6 5 ’ 19 4.84 5 - 15 4.83 

n 0,5 

“ 

■I Workday Weekend 

■ Ave rage Journey Time ■ Ave rage Runni ng Time 
■ Average Del ay time 



Fig. 3 Travel time & delay of off-peak hour (N-S) trip. 


Travel time & delay of Off-Peak Hour 
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■ Average Journey Time ■ Ave rage Running Time 

■ Average Delay time 

Fig. 4 Travel time & delay of off-peak hour (S-N) trip. 

In workdays for the South to North (S-N) trip the average 
delay time was the maximum (0.56 minutes) during peak hours. 
And during the off-peak hours of weekend days, the delay time 
was the minimum (0.30 minutes) for the same South to North (S- 
N) trip. Even the average delay time (0.40 minutes) during the 
peak hours in the weekend days (Fig. 1) is greater than the 
average delay time during the off-peak hours of working days 
(Fig. 3). 

Normally in peak hours, the number of vehicles remains 
higher than the off-peak hours. So that means the number of 
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vehicles can be the prominent factor for the increased amount of 
delay time during the peak hours in the workdays. If the delay 
time of the study area becomes perfectly equivalent to vehicle 
flow, then it will mean that the total number of vehicles is the 
number one factor for delay time in this road. And it will also 
mean that if the number of vehicles increases every time the 
delay time also increases in this road too if the other factors 
remain the same. 

4.2 Vehicle flow analysis and related travel time-delay time 
comparison 

To analyze the traffic flow behavior of the study area the 
vehicle flow chart was prepared. It also explains the possible 
cause of journey time and delay time differences during peak and 
off-peak hours for both workdays and weekend days. The 
comparative vehicle flow for both directions of the road during 
the peak and off-peak hours in workday and weekend is shown 
in Fig. 5 and Fig. 6. 
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Fig. 5 Vehicle flow of N-S trip. 
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Fig. 6 Vehicle flow of S-N trip. 

Workday had the highest Vehicle flow in peak hour for both 
directions of a vehicle trip. Here, the S-N trip shows a bigger 
number of vehicle flow in this study area. Most importantly, the 
total number of vehicles during peak hours of workdays was 833, 
and 913 which are greater than the off-peak hour’s vehicle flow 
in workdays. These relations remain the same for both North to 
South (N-S) and South to North (S-N) trips. That means in this 
study area the total number of vehicles is always higher during 
the peak hours and hence the traffic pressure also remains higher 
during the peak hours of both workdays and weekend days and 
for both side of the trip. 

The journey time and delay time data also show the exact 
fact of increased traffic pressure during the peak hour of both 


workdays and weekends in both directions of the roads (Fig. 1- 
4). So in the study area, there is a direct relation of vehicle flow 
number and the average delay time. According to the data shown 
in Fig. 1-6, as the total number of vehicles increases in peak 
hours, the delay time also increases as result. 

4.3 Spot speed analysis 

Spot speed analysis was performed to determine the speed 
distribution of different vehicles (Easy bike, Mahindra, 
Motorcycle, Bus, Truck) available in the traffic stream of this 
study area. Spot speed frequency distribution curve and 
cumulative distribution curve have been analyzed for both 
intersections (Religate and Mohosin More) to find the vehicle 
types which are directly responsible for the increased delay time 
in this study area. 

4.3.1 For Relegate intersection 

The collected speed datasets were grouped into several 
speed-class interval sections. Within these sections, the number 
of vehicles was counted to prepare the speed frequency 
distribution table. From this table, the arithmetical mean speed 
and a histogram of speed-class intervals against the percentage 
frequency of vehicle numbers were developed (Fig. 7). Here the 
speed represents the average speed of all the collected vehicle’s 



Fig. 7 Spot speed frequency distribution curve for Religate 
intersection area. (Source: Field Survey) 



Fig. 8 Cumulative speed distribution curve for Religate 
intersection road area. (Source: Field Survey) 

The average speed of vehicles in the Religate intersection 
area was 26.0 K.P.H. with a standard deviation of 12.66 which 
means speed fluctuates a lot during traffic congestion. From Fig. 
7 it is clear that the speed of vehicles does not tend to cluster to 
the mean speed. But the frequency increases if the speed deviates 
slightly from the mean speed value. This indicates that there are 
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It means half of the vehicles travel faster than this 50 th 
percentile speed and half of the vehicles travel slower 
than this. 

iii. The 85 th percentile for the Religate intersection is 41.20 
K.P.H. which represents that under a free-flowing 
condition 85% vehicles go in this speed or slower than 
this. And only 15% of vehicles travel faster than this 
speed. The 85 th percentile speed is commonly used in 
traffic engineering and traffic safety guidelines as it is 
considered as the quickest safe speed for any traffic in 
a road under a free-flowing condition [19]. 

iv. The 98 th percentile speed for the Religate intersection 
area is 47.50 K.P.H. That means 98% of the vehicles 
travel equal or less than this speed. This speed is 
normally considered as the design speed of any road 
[ 21 ]. 

Table 3 reveals that Mahindra and Easy bike have the 
highest coefficient of variation with a maximum speed of 30 
K.P.H and 23 K.P.H. respectively. It means that they could move 
faster if they don’t face any congestion. Mahindra and Easy bike 
account for a very large portion of traffic composition on the 
road. The minimum speed of Mahindra and Easy bike go lower 
than the 15 th percentile speed for this section. 


Table 3 Spot Speed statistics of different types of vehicles on Religate intersection. 

Speed of Vehicles 

v cmcic type 

Mean 

Maximum 

Minimum 

Standard 

Deviation 

Coefficient of 

variations 

Easy bike 

12.31 

23.00 

4.00 

5.95 

.48 

Mahindra 

17.40 

30.00 

5.75 

8.90 

.51 

Motorcycle 

33.33 

53.20 

12.41 

12.41 

.37 

Private Car 

34.35 

47.50 

15.20 

10.21 

.30 

Bus 

32.59 

45.60 

13.50 

9.63 

.30 

Truck 

26.31 

43.29 

12.50 

10.36 

.39 


a significant amount of vehicles travel both at high speed and low 
speed in this section. 

4.3.1.1 Cumulative distribution curve 

The cumulative speed distribution curve for the relegate 
intersection area has been prepared from the speed-class interval 
and frequency table. To prepare the curve the speed-class 
intervals are plotted in the X-axis the cumulative percentage of 
speeds for each class intervals are plotted on the Y-axis. Then 
from the curve, the related percentile speed is identified (Fig. 8). 
The related percentile speed for Fig. 8 are as follows: 

i. The 15 th percentile speed for Religate intersection road 
is 10 K.P.H. That means 85% of vehicles tends to go 
faster than 10 K.P.H. in this section. And hence vehicles 
with less than or equal 10 K.P.H. are responsible for 
increasing the delay time and congestion significantly 
in this area. So the 10 K.P.H speed limit should be 
prohibited in this section. This will reduce congestion 
and the unnecessary overtaking tendency of the other 
drivers in this section. 

ii. The 50 th percentile speed for this section is 23.6 K.P.H. 
which the median speed for the observed data set. It 
represents the average speed of the total traffic stream. 


It means that if queuing is formed during congestion, this 
traffic mode makes the congestion more acute as other vehicles 
need to overtake opportunities while approaching to move. The 
direct field survey data (Fig. 9) also upholds the results of Table 
3. The minimum spot speed of Easy and Mahindra found lower 
than the 15 th percentile speed in this intersection with a speed of 


only 4 K.P.H and 5.75 K.P.H respectively. According to the field 
survey data, Easy bike and Mahindra often stops in this 
intersection to collect passenger which creates an unwanted 
queue. This consequently results in traffic congestion and travel 
time delay. 
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(e) Mean Speed of Motorcycle (f) Mean Speed of Private Car 


■ 12.31 km/h 

■ 17.40 km/h 

■ 26.31 km/h 

■ 32.59 km/h 

■ 33.33 km/h 

■ 34.35 km/h 


Fig. 9 Mean speed map for different vehicles in the Religate More intersection. 


4.3.2 For Mohosin More intersection 

The mean speed of Mohosin More intersection was 25.7 
K.P.H. with a standard deviation of 12.3. That means most of the 
speed varies from 13.4 K.P.H to 38 K.P.H. The spot speed 
frequency distribution curve shows that the speed of vehicles was 
not also showing cluster tendency to the mean speed in the 
Mohosin More intersection. But the maximum vehicle speed was 
very close to the mean value. However, the frequency decreases 
if the speed deviates from the mean speed value. 

The curve indicates that there were a significant amount of 
vehicles that traveled near the mean speed value though the 


overall number of vehicles with high speed and low speed were 
also significant in this intersection. 

4.3.2.1 Cumulative distribution curve 

Using the Spot speed frequency distribution curve’s data the 
cumulative speed distribution curve has been prepared. From this 
curve, the 15 th , 50 th , 85 th ’ and 98 th percentile speed has been 
determined that helps to determine those which travel speed is 
lower than the desired speed and hence they can create extra 
traffic pressure and delay time as a consequence. 
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Mean Speech 25.7 K.P.H. 



Fig. 10 Spot speed frequency distribution curve for 
Mohosin More road area. 



Fig. 11 Cumulative speed distribution curve for Religate 
intersection road area. 


The related percentile speeds for Mohosin More intersection 
area are as follows: 

i. The 15 th percentile speed for the Mohosin More road 
intersection is also 10 K.P.H. which is just similar to the 
Relegate intersection area. That means 15% of the total 
vehicles travel equal to or lower than this speed. This 
speed limit should be prohibited in this section too as 
this speed limit increases delay time, congestion time, 
and overtaking tendency to the other drivers. 

ii. The 50 th percentile speed for this section is 24.5 K.P.H. 
That means 50% of vehicles travel faster than this speed 
and the rest of the 50% travels slower than this speed 
limit. 

iii. The 85 th percentile speed for this road section is 41.38 
K.P.H. That means in a free-flowing condition 85% of 
all vehicle travel within this speed limit. Normally this 
speed limit is considered as the maximum allowable 
safe speed for any road. 

iv. The 98 th percentile speed for Mohosin More Road area 
is 48.75 K.P.H. This speed can be considered as the 
design speed for this road. 


Table 4 reveals that Mahindra and Easy bike have the 
highest coefficient of variation with a maximum speed of 29.50 
K.P.H and 24.29 K.P.H. respectively. It means that they could 
move faster if they don’t face any congestion. Mahindra and 
Easy bike account for a very large portion of traffic 
composition on the road. The minimum speed of Mahindra and 
Easy bike go much lower than the 15th percentile speed for this 
section. It means that if queuing is formed during congestion, 
this traffic mode makes the congestion more acute as other 
vehicles need to overtake opportunities while approaching to 
move. 

The direct field survey data (Fig. 12) also approves the 
results of Table 4. Just like the Religate intersection, in the 
Mohosin More intersection, the minimum spot speed of Easy 
bike and Mahindra found lower than the 15 th percentile with a 
speed of only 5 K.P.H and 4 K.P.H respectively. According to 
the field survey data, Easy bike and Mahindra tends to stop in 
these intersections to collect passenger which creates unwanted 
queue in this study area. This consequently adds extra traffic 
pressure and travel time delay as result. 


Table 4 Spot speed statistics of different types of vehicles on Mohosin More intersection. 

Speed of Vehicles 

vcmuic iypc 

Mean 

Maximum 

Minimum 

Standard 

Deviation 

Coefficient of 

variations 

Easy bike 

13.32 

24.20 

5.00 

6.09 

0.46 

Mahindra 

17.45 

29.50 

4.00 

9.13 

0.52 

Motorcycle 

31.31 

53.25 

15.50 

12.81 

0.41 

Private Car 

34.77 

48.75 

16.20 

10.63 

0.31 

Bus 

30.57 

43.20 

12.25 

9.24 

0.30 

Truck 

26.76 

42.75 

12.50 

10.16 

0.38 
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(e) Mean Speed of Motorcycle (f) Mean Speed of Private Car 


■ 13.32 km/h 

■ 

17.45 km/h 

■ 

26.76 km/h 

■ 30.57 km/h 

■ 

31.31 km/h 

■ 

34.77 km/h 


Fig. 12 Mean speed map for different vehicles in the Mohosin More intersection. 


5. Conclusion 

A comparison among travel time, running time, delay in 
peak, and off-peak hours on different days of a week is out in this 
study. From the result and analysis conclusion has been made 
that workdays have high travel time & delay compared to the 
weekends. The study shows the relation between traffic flow and 
delay time, as the traffic flow increases the delay time also 
increases. On the other hand, the reasons behind the delay time 
are also found out in Religate and Mohosin More intersections. 
Vehicles with less than or equal 10 K.P.H. are responsible for 
increasing the delay time and congestion significantly in these 
intersections, as it is the 15 th percentile speed for these 
intersections. Easy bike and Mahindra have a minimum speed of 
5 K.P.H and 4 K.P.H respectively at Mohosin More intersection 
and for Religate intersection the minimum spot speed of Easy 
and Mahindra are 4 K.P.H and 5.75 K.P.H respectively. As the 


minimum speed of Mahindra and Easy bike are much lower than 
the 15 th percentile speed for these sections, they are responsible 
for increasing the delay time and congestion significantly on this 
road. 
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ABSTRACT 


A quest for physics-based understanding of engineered metamaterials inspired numerous researchers to extract intricate features 
such as guiding and filtering elastic waves, wave focusing, topological insulation, cloaking etc. A traditional metamaterial is composed 
of a heavy core along with some other materials of dissimilar mechanical properties (e.g. stiffness, density). It is well-established that 
the frequency band gap can be formed by introducing the desired inconsistency of material properties between stiffer resonators and 
adjacent embedding matrices. Frequency band gaps are fundamental requirements of many engineering applications such as vibration 
control, noise mitigation, and energy concentrations. Hence, advanced researches are being carried out continuously to understand the 
control parameters (e.g. bandwidth, starting and ending frequencies) of the frequency band gap. In this article, a mass-in-mass 
metamaterial using elliptical anisotropic resonators are considered to investigate the influence of resonators’ geometric factors and 
volume fraction on band gap parameters. While the elliptical resonator is splitted diagonally in one or both opposite ends to analyze 
the influence of volume fraction of the resonators on stop bands, a second set of half-circular resonators are also investigated to analyze 
the impact of resonator parameters on frequency bands. 


Keywords: Acoustic metamaterial; Band gap; Split-ring resonator; Volume fraction; Wave filtration. 
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1. Introduction 

From the last couple of decades, electromagnetic 
metamaterials [1],[2] are in concentration to phononic 
researchers to explore photonic band gaps, a range of frequency 
where electromagnetic waves cannot propagate [3]. The concept 
of electrodynamics metamaterial can be used in designing 
elastrodynamic metamaterials since the physical understanding 
of electromagnetic wave propagation in photonic crystals is quite 
similar to the stress (elastic) wave propagation in phononic 
crystals [4],[5]. Mass-in-mass systems are frequently proposed 
under elastodynamic problems to predicatively manipulate the 
frequency band gaps in metamaterials with engineered volume 
inclusions [6],[7]. Creating local resonance vibration modes is 
one of the important features of such acoustic metamaterials. 
Applications of such acoustic metamaterials have been 
envisioned for waveguiding, vibration control, sound isolation, 
sensing, energy harvesting etc. [8]-[12]. 

Band gap and band structure manipulation are increasingly 
important in acoustic metamaterials to branch out its 
applicability. It is well-established that in any metamaterial 
frequency band gap are results from either local resonance or 
Bragg scattering. Fow-frequency sound can be controlled by 
introducing locally resonant components into phononic crystals 
[13]-[17], whereas conventionally high-frequency stop bands 
can be formed by multiple scattering (Bragg) of the periodic 
inclusion [18],[19]. In recent years, researchers are attempting to 
manipulate the ability of the metamaterial by artificially 
designing the system. It has been proven that the number of stop 
bands can be increased by adding additional local resonators in 
the system [20],[21]. Most recently a multi-scale mass-in-mass 
model (MMM) with split-ring resonators [22] is proposed not 


only to obtain a wider frequency band but also to obtain multiple 
band gaps in both sonic and ultrasonic frequency regions. The 
concept of using split rings in a metamaterial system has 
primarily emerged in photonics research. Movchan et al. [23] 
used split-ring resonators to control electromagnetic bands in 
two-dimensional photonic structures. Many other researchers in 
photonic domain found split ring useful in manipulating 
electromagnetic waves for specific purposes. 

Inspired from these researches, a symmetric split-ring 
metamaterial system with multiple resonators is considered 
herein as the model offers a high level of anisotropy. Note that, 
anisotropy in a unit cell is one of the key factors of dispersive 
behavior and band gap formation in a mass-in-mass metamaterial 
system. Although the size effect of the resonators on the 
respective band gaps was reported earlier [24], the geometric 
(shape)/volumetric dependency of the resonators on frequency 
bands is still to explore. Hence, the primary objective of this 
study is to investigate the possibility to systematically 
manipulate the frequency band gap by the alteration of the 
resonator’s geometry or volume fraction. In addition, two 
different sets of resonators in a unit cell are also investigated to 
understand the influence of resonator setting on the frequency 
band gap. 

2. Model and Computational Approach 

2.1 Model Configurations 

Several two-dimensional structures with multiple resonator 
systems (Ref. Fig. 1) are designed to understand the 
geometric/volumetric dependency of the local resonators on 
frequency bands. Initially, the unit cell is composed of a heavy 
core (Rl) of diameter 0.1414 inches embedded in a circular ring 
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(R2) with an outside diameter of 0.2828 inches. A softer material 
(Ml), say, Rubber, is used to seal the space between the heavy 
core and the circular ring. A similar mass-in-mass system was 
also proposed by Huang et al. [21]. However, an elliptical ring 
(R3) is placed at the center of the unit cell to further enhance the 
anisotropic property of the system (Fig. 1 (a)). The elliptical ring 
is rotated at 450 about the center of the unit cell to allow identical 


wave pathway along both x- and y-directions. A much stiffer 
epoxy material M2 compared to Ml is considered to fill in the 
elliptical free space. In order to allow continuous mismatch in 
material stiffness between components in the system, a polymer 
material (M3) is considered to complete the outer boundary of 
the unit cell. The outside dimension of the unit cell is measured 
as one (1) square inch. 






Fig. 1 Projected models to analyze the geometric influence of local resonators on frequency bands. Models (la-lc) with 
resonators R1-R3 are considered as set-1 and rest (Id-If) as set-2 with R1-R4 resonators. 


While designing, the elliptical resonator is splitted about 
0.115 inches along its major axis in both one (Fig. 1(b)) and two 
(Fig. 1(c)) directions, respectively. Such splitting of the 
resonators not only decreases the volume fraction of the 
resonators in the unit cell but it also allows additional degrees of 
freedom in the system, which tends to alter the frequency bands. 
The splitting of the resonators was made at the ends of the major 
axis of the elliptical geometry to ensure similar wave propagation 
path along x- and y- directions. In order to analyze the geometric 


influence of the resonators more rigorously, especially the effect 
of resonator setting, a pair of half-circular rings (R4) are added 
symmetrically (Figs. 1(d)-1(e)) in all three models described 
previously (Figs. 1(a)-1(c)). Unit cells with resonators R1-R3 
are considered as resonator set-1 (Models la-lc); whereas set-2 
is referred to as models with resonators R1-R4 (Models ld-lf). 
In this study, the thickness of all the rings is considered as 0.037 
inches. Material properties and dimensions of the unit cell 
components are listed in Table 1. 


Table 1 Properties of the components of the unit cell. 


Component 

Name 

Outer Dimensions 
(inch) 

Stiffness 

(Pa) 

Density 

(kg/m 3 ) 

Poisson’s 

Ratio 

Ml 

Diameter-0.2121 

10 e6 

980 

0.49 

M2 

Major Radius - 0.4 
Minor Radius - 0.2 

2.5e9 

1250 

0.38 

M3 

1X1 Square 

0.5e9 

1050 

0.49 

R1 

Diameter - 0.1414 

13e9 

11310 

0.435 

R2, R4 

Diameter - 0.2828 

100e9 

2950 

0.31 

R3 

Major Radius - 0.435 
Minor Radius - 0.235 

100e9 

2950 

0.31 
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2.2 Numerical Simulation 


The mathematical formulation of infinitely repeating split 
ring unit cell proposed in this study requires the implementation 
of Bloch theorem. Since it is challenging to perform analytical 
derivation of such a material system, the widely accepted Finite 
Element Method has been utilized to solve this problem. While 
calculating the dispersion relation, the complete structure is 
considered infinite in both x- and y- directions by arranging the 
unit cells or the Representative Volume Elements (RVEs) 
periodically. To implement this concept, the Bloch-Floquet 
periodic boundary condition is applied at all boundaries of the 
RVE. These boundary conditions are based on the Floquet theory 
that can be applied to the problem of vibrations with small- 
amplitude in spatially periodic structures. The theory states that 
the solution can be sought in the form of a product of two 
functions. One follows the periodicity of the structure, while the 
other one follows the periodicity of the excitation. The problem 
can be solved on a unit cell of periodicity by applying the 
corresponding periodicity conditions to each of the two 
components in the product. Generalized wave equation in the 
composite material can be written as 

Cijkl(Xm)[u k> i(x m , t) + U l k (x m , t)] + fi(x m ) 

P (Xm)iii(x m i t) 


where the constitutive matrix containing material properties 
and the density of the system are the functions of space (x lt x 2 ). 
i,j, k,l & m takes values 1, 2 & 3. Standard index notation is 
used throughout this manuscript. Let the body force /(x m )be 
constant. 

The unit cells are repeated in both directions and the solution 
can be assumed in terms of Bloch-Floquet solution as discussed 
in the previous paragraph. Assuming there is no periodicity along 
x 3 direction & decoupling the phase component, the 
displacement solution can be viz. 


Ui(x m ,t)= 'Y^A i nln2 exp(ik m x m ). 


( 2 ) 


exp(iG m x m ). exp{ik 3 x 3 ) . exp(-icot) 

Where, k m is wave number along 772-th direction and G m is 
a component of the reciprocal lattice vector along m-th direction. 
Here, m takes values 1 & 2. G m can be expressed as G m = 
2nn m /D m , where, D m is the periodicity of the cells in m-th 
direction. The A l nln2 is the amplitude of the wave modes for 
particle displacement along i and nn 2 are the integer numbers 
between — oo to + oo. After substituting Eq. (2) in Eq. (1) the 
Bloch equation with the Bloch operator can be obtained as 
follows 


0) 


! POm) ^ ^ A l nln2 exp (i(k m + G m )x n 


712 7ll 


2 

+ G m ) 2 S m ^j ^ ^A k n i n2 exp(i(ki + G t )x{) 
+ / / A l n i „ 2 exp(i(k k + G k )x k ) = 0 


(3) 


The above equation is a Bloch eigenvalue problem. The Eq. 
(3) is then multiplied with Bloch operator with Bloch 
transformed weighting factor and integrated over the whole 
domain of the body and the equation was transformed to its weak 
form. Periodic boundary conditions are applied around the unit 
cell and weak form of Bloch equation is solved only within the 
irreducible Brillouin Zone. Further number of amplitude in the 
Eq. (3) is reduced for each wavenumber (k) point. Thus the 
n! & 72 2 are reduced from infinity and the truncated set of Bloch 
mode expansions were used in the solution method. A suitable 
choice of reduced-order basis was made based on the high 
symmetry points that characterize the periodic lattice. Next, the 
Finite Element discretization was performed. Triangular 
isoparametric elements were used in the simulation. Element 
sizes were determined based on a series of convergence study 
and the minimum wavelength occurred in the material. The sizes 
of the elements were kept to a minimum of 1/10 of the 
corresponding minimum wavelength that occurred in any 
material type, respectively. The Bloch displacement amplitudes 
were discretized using isoparametric shape function (N t (x)) 
suitable for triangular elements for each combination of ^ & n 2 
in their truncated series as follows 

3 

A n in2=Y J N iW A i (4) 

i =1 

Applying the discretization equations and periodic boundary 
conditions the weak form of Bloch equation reduces to an 
algebraic eigenvalue problem[K(/e) — co 2 M]V = 0, where V is 
the discrete Bloch amplitude vector which is periodic within the 
unit cell. The K(k) and M are the global stiffness and mass 
matrices, respectively obtained by integrating the element level 
matrices in proper order. Detail expressions for K and M can be 
found in reference. The solution of the eigenvalue problem 
provided the dispersion curves for the proposed periodic media. 

2.3 Validation of Computational Approach 

The studied model configurations in this work are quite 
complex. Hence, it is extremely challenging to perform the 
mathematical formulation and analytically solve the proposed 
models. Hence, a simple mass-in-mass system (Ref. Fig. 2) is 
studied herein using the Finite Element tool to validate the 
solution methodology. 


90 D 
33.9 D 
25.3 D 



105 


Fig. 2 Single-resonator mass-in-mass system with square 
periodicity (dimensions are in mm). 
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In this simple system, a unit cell is designed which is 
composed of four components. A significantly stiff material, lead 
(component 1), is chosen as a circular resonator [4],[25]. To 
avoid evanescent wave modes, the core resonator is coated with 
softer rubber-like material (component 2). The rubber-lead 
arrangement is then embedded into a relatively stiffer component 
(epoxy). The lead-rubber-epoxy combination is synonymous 
with a well-known mass-in-mass engineered material system 
proposed by the earlier researchers [20]. It is mathematically 
proven that such an engineered system provides frequency band 
gaps (stop bands) by virtually creating negative bulk modulus 
and negative mass density in the structure [21],[26]. In 
engineered material systems, the negative responses result from 
the mismatch of material properties between adjacent 
components. Hence, to induce further property mismatch, the 
mass-in-mass system is placed inside a square unit element and 
the created gap is filled with a relatively softer material 
(polyethylene plastic). 

The complete lead-rubber-epoxy-polyethylene design is 
considered as the Representative Volume Element (RVE) in this 
section. The unit cell is a 105 mm square and diameters of the 
lead, rubber, and epoxy components are considered as 25.4 mm, 
33.9 mm, and 90 mm, respectively. The material property of the 
unit cell components is listed in Table 2. 


Table 2 Properties of the components enclosed in the unit cell. 



Lead 

Rubber 

Epoxy 

Polyethylene 

Young’s 
Modulus (Pa) 

13e9 

10e6 

3.5e9 

0.7e9 

Density 

(kg/m 3 ) 

Poisson’s 

Ratio 

11310 

980 

1250 

1050 

0.44 

0.49 

0.38 

0.49 



ka/n 


(a) (b) 

Fig. 3 (a) Dispersion relation (b) Mass density plot for the 
single-resonator mass-in-mass system. 

It is well known that the frequency band gaps are the result 
of negative mass density and induction of resonance in local 
microstructures of metamaterials [20] , [21 ] , [ 14] , [24] . To validate 
the above statement and adopted FEM solution methodology, a 
multi-layered mass-in-mass model (Ref. Fig. 2) is considered. 
This study is carried out for a frequency range of 0 to 10 KHz. 
For illustration purposes, the response is reported between 0 to 4 
KHz in Fig. 3. Computed effective mass densities are normalized 
by the density of the stiffest component (lead) of the cell. The 
effective mass density of the unit cell is calculated by Eq. (5) 
using long wavelength assumption. 


_ / Gapdr 

Pap fUapdr 


(5) 


where a, /? = 1,2 and a a p and ii a p are the local stress and 
acceleration quantities, r denotes the external boundary of the 
unit cell. 


The dispersion curve for the single resonator metamaterial 
within the first Brillouin zone is shown in Fig. 3(a). Using the 
definition of mass density in Eq. (1), it has been found that the 
mass density of the system stays negative between the frequency 
range of 2.27 - 2.71 KHz (Ref. Fig. 3(b)). Alternatively, a stop 
band is observed (Ref. Fig. 3(a)) at a frequency range between 
2.28 - 2.74 kHz (considering T-X directional waves). The same 
band gap continues for the entire T-X-M-T directional waves. 
Since the mass density of the metamaterial (calculated using Eq. 
(1)) is found negative within the frequency range of the stop 
band, it can be concluded that the analytical and numerical 
outcomes are in close agreement. This phenomenon is well 
reported in the literature [20],[21],[14],[24],[27],[28]. 

3. Results and Discussion 

In order to investigate the influence of resonator geometry, 
six (6) models are proposed as described in Fig. 1. It is 
hypothesized that resonator geometry as well as the volume 
fraction of resonators in the unit cell and type of resonator 
arrangement may have a considerable impact on frequency 
bands. The dispersion relation of the irreducible Brillouin zone 
(see Fig. 4(b)) for model la (Ref. Fig. 1(a)) is shown in Fig. 4(a). 
Note that, the wave vector is normalized by the length of the unit 
cell ‘a’. 
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(b) 

Fig. 4 (a) Dispersion curve and band gap representation 
for model la. (b) First Brillouin Zone. 

In Fig. 4(a), almost all the frequency bands seem highly 
dispersive, whereas a strong weakly-dispersive (straight) region 
is observed in the first frequency band, and a stop band is noticed 


between 7 th and 8 th frequency bands. The primary focus in this 
section is to study the transformation of the band gaps due to the 
alteration of resonator geometry and settings, For comparison, 
the magnified region of the band gaps for all the models are 
plotted in Fig. 5 (set 1) and 6 (set 2). 




Normalized wavevector (ka/n) 


Fig. 6 Magnified dispersion curve representations for resonator set-2 (models Id-If). 


A band gap with a bandwidth of 562 Hz is observed between 
21.598 KHz to 22.16 KHz for the unit cell with a center ball, 
circular ring, and closed elliptical ring resonators (Fig. 5(a)) and 
Table 3). However, the bandwidth is reduced to 402 Hz (Fig. 
5(b)) and 241 Hz (Fig. 3(c) and Table 3) through the splitting of 
the elliptical resonator in one end and both ends, respectively. A 
similar trend is also noticed for other sets of models where two 
additional half-circular resonators are placed symmetrically. The 
longest stop band is recorded for model Id with a bandwidth of 


792 Hz and reduced to 661 Hz and 528 Hz, for models le and If, 
respectively. 

Volume fractions of the resonators in the respective unit cell 
are listed in Table 3. It is evident that for each set of models, the 
stop bandwidth can be increased/decreased with the 
increment/decrement of the volume fraction of the resonators, 
and it follows a linear pattern. Such a linear pattern is true for any 
particular set of resonators and can be shifted significantly 
through the alteration of the resonator set as found in Fig. 7. This 
finding implies that the same stop band can be maintained with 


89 































































R. Ahmed, H. Ahmed and S. Banerjee/JEA VoL 01(03) 2020, pp 85-93 


different volume fractions of the resonators and there is also a fraction. In both cases, the selection of resonators configuration 

possibility to obtain longer band gaps with unchanged volume in the unit cell would be the key challenge. 


Table 3 Obtained band gaps and volume fraction of resonators in the unit cell for studied models. 


Model # 

Band Start 

Band End 

Band Gap 

Volume Fraction 

(Hz) 

(Hz) 

(Hz) 

(%) 

a 

21598 

22160 

562 

11.38 

b 

21608 

22010 

402 

11.17 

c 

21623 

21864 

241 

10.96 

d 

21454 

22246 

792 

14.13 

e 

21453 

22114 

661 

13.91 

f 

21457 

21985 

528 

13.70 



Fig. 7 Stop band manipulation through the alteration of 
resonators volume fractions and setting in unit cell. Each bullet 
defines the obtained stop bandwidth with respect to the volume 
fraction of the resonators in the corresponding proposed model 
(la-lf, see Fig. 1). 

Fig. 7 shows the influence of the volume fraction of the 
resonators on OVERALL stop bandwidth. Here, OVERALL 
means, these stop bands can be found for wave propagation 
through the model along any direction. While the overall band 
gap of a structure is important for filtering wave propagation in 
a three-dimensional space, guided wave propagation (along a 
specific direction) in a structure is not uncommon. Fig. 8(a) 
represents the stop band widths for both set-1 and set-2 if the 
wave enters the structure as a guided x-directional (r-X) wave. 
Fig. 8(b) represents a similar phenomenon in the case of the y- 
directional (X-M) wave. It can be seen that the stop bands are 
heavily influenced by the volume fraction of the resonator for 
wave propagation in T -X direction, which is similar to the overall 
response reported in Fig. 7 However, in the case of wave 
propagation in X-M direction, the effect of volume fraction isn’t 
much significant, where a little increase in stop bandwidth is 
identified (Fig. 8(b)). 

In a bandgap region, it can be seen that two bands from the 
dispersion curve take a role in forming the stop band (Ref. Fig. 
4(a)). To understand the influence of resonator volume fraction 
on those band gap forming frequency bands, a further 
investigation is performed and the results are presented in Figs. 
9-11. Fig. 9 presents the shift in the lower and upper band of the 
overall bandgap for different volume fractions. Figs. 10 and 11 
present similar outputs for x-directional and y-directional guided 


waves, respectively. In all these figures, it can be noticed that the 
volume fraction of the resonators in the unit cell possesses almost 
no influence on the lower band (start frequency) of the stop band 
for both sets (a-c and d-f) of models. 




Fig. 8 Influence of volume fraction on stop bandwidth for 
wave propagation along (a) T-X or x-direction, and (b) X-M or 
y-direction. 

However, the upper bands of the bandgap show a significant 
dependency on the volume fraction. Note that, for X-M 
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directional wave incidence (Fig. 11) both the upper and lower 
band of the band gap show very little reliance on volume fraction. 
Henceforth, the bandgap width also remains almost unchanged. 
Upon investigation of the mode shape, it has been found that the 
frequency bands forming the stop bands along X-M direction are 
predominantly the result of Bragg scattering, where deformation 
of the matrix element in the unit cell controls the formation of the 
band structure. In such instances deformation of the ring 
resonators is quite negligible. Hence, a change in the volume 

Models a-c 



fraction of the ring resonators making a negligible effect on the 
band structure or band gap along X-M direction. On the contrary, 
the upper band to form the frequency band gap along Y-X 
direction is the result of the combined effect of Bragg scattering 
and local resonance of the metal resonators. Such observation is 
also true for the upper bands of the overall stop bands in Fig. 9. 
Hence, a change in the volume fraction of the resonators shifted 
the upper bands for both overall and Y-X directional waves. 


Models d-f 



Fig. 9 Shift in start and end of the overall band gap due to the change in volume fraction 

Models a-c Models d-f 



Volume fraction (%) 

(a) 



Fig. 10 Shift in start and end of the band gap along Y-X due to the change in volume fraction 

Models a-c Models d-f 




Fig. 11 Shift in start and end of the bandgap along X-M due to the change in volume fraction 
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In Figs. 9-11, it can be seen the influence of the resonator’s 
volume fraction on stop bandwidth. However, because of the 
different scaling in both x- and y-axis it can be a little difficult to 
identify how much the volume fraction is taking part in 
controlling the bandgap width for various situations. For 
example, in Figs. 9-10, all the presentations look almost the 
same. In that context, the volume fraction possesses almost equal 
influence for both the set of models and for both overall and x- 
directional waves. To have a clearer understanding, a percentage 
increase in bandgap width due to the increase of volume fraction 
is plotted in Fig. 12. It can be seen that using the model set-1, 
maximum increase (about 135%) can be obtained in the case of 
overall bandwidth (ref. Fig. 12(a)). Whereas, almost negligible 
increase is visible in X-M directional wave. A similar trend can 
be visible for model set-2 as well (ref. Fig. 12(b)). However, the 
maximum increase for set-2 in the overall bandgap is about 50%, 
which is significantly lower than the gain reported for set-1 in 
Fig. 12(a). 


containing split-ring resonators are investigated to understand 
the influence of resonator setting and volume fraction on the 
frequency stop bands. It has been found that, in any resonating 
setup, width of the frequency bandgap can be linearly 
manipulated through the alteration of volume fraction of 
resonators in the unit cell. However, the linear pattern can be 
shifted sidewise by changing the resonator set by introducing or 
taking off a resonator from the unit cell. It has been noticed that 
lower band of the stop band remains unaffected with the change 
of volume fraction and the upper band can be shifted 
upward/downward through the increase/decrease of volume 
fraction of the resonators. This study also confirms that, with a 
lower number of resonators it is possible to manipulate the 
frequency band gap heavily by changing the volume fraction, in 
comparison to a higher number of resonators. In totality, this 
study concludes that, in a metamaterial, frequency bands are 
strongly sensitive to resonator setting and volume fraction in a 
unit cell. 
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Fig. 12 Percentage increase in band gap width for models 
(a) la-lc (set 1), (b) ld-lf (set 2) 



4. Conclusion 

Manipulating frequency band gaps by designing novel 
configurations of metamaterials is one of the key areas of interest 
in recent days. In this purpose, mass-in-mass and layer-in-layer 
models were proposed by the previous researchers. However, 
dependency of frequency band gaps on resonator setting and 
volume fraction of the resonators in a unit cell were still to 
investigate. In this work, multi-scale mass-in-mass systems 
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ABSTRACT 


Power requirements are growing day by day, and more power plants are being constructed all over the world. Now the goal has 
been to look for more sustainable sources of electricity. Sustainability of power plants is a complicated concept and depends on various 
criteria and sub-criteria. By evaluating them separately creates a complex problem. For this here Multi-Criteria Decision Analysis 
(MCDA) methods are used for overall assessment and make a sustainability index for seven mostly used power plants. Consideration 
is taken to both renewable and non-renewable sources. The goal is divided into three basic criteria (i.e. technology, safety & 
sustainability, economy) and each criterion is further sub-divided into different sub-criteria. The data is collected from various sources 
and then analyzed using AHP and PROMETHEE methods. The result indicate renewable sources are typically advantageous over non¬ 
renewable sources. In certain cases, nuclear has some benefits over other non-renewable energy sources. 


Keywords: Sustainability; Power plants; MCDA; AHP; PROMETHEE. 
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1. Introduction 

Sustainability can be termed as a process that satisfies the 
present needs of people without hampering the generation 
process for future generations. The sustainability concept is 
based on three main pillars: social, economic and environmental. 
The term sustainability was first brought up by the UNESCO in 
1987 as a part of boosting the development process. There are 
goals and visions for this sustainability which take the society to 
a brighter future [ 1 ] . 

As the need of electrical power is day by day increasing 
there should be emphasis given on the power sector. Power plant 
building without any pre-analysis or pre-planning is pretty much 
unwise by considering future needs of power. As power can be 
generated from different sources it is wise to analyze which 
power source to be chosen for more effective performance of the 
power plant not just in the present but also in the future. So as for 
this, the sustainability term comes in work. The sustainability of 
power plants means selecting the appropriate power plant 
considering different condition and selecting a power plant for 
both present and future needs. 

In Bangladesh defining and understanding the term 
sustainability for power plants has become very important 
nowadays because a rumor has been built up whether the 
construction of new power plants like Matarbari coal-based 
power plant, Rampal coal-based power plant, Rooppur nuclear 
power plant etc. are sustainable or not. And why these fuel 
sources are being chosen instead of other sources likewise gas, 
oil, wind, solar, etc. The sustainability index can answer this sort 
of question about their importance in the present and future 
periods. As being a developing country Bangladesh should keep 
eye on the sustainable development of power source as the 
requirement of power will increase exponentially day by day. 

For developing a sustainability index for different types of 
power plants at first the goal is to define sustainability in terms 
of parameters which can be measured and expressed by data. To 


do that sustainability is divided into three main criteria as 
technology, safety and sustainability, and economy [2]. These 
three criteria are then subdivided into eleven sub-criteria. Now to 
make a sustainability index, multi-criteria decision analysis 
(MCDA) method is used. MCDA is the collection of decision¬ 
making methods that helps a decision-maker to select a decision 
based on different data. The data are collected from trusted 
sources and after analyzing the data a goal can be achieved. 
Different types of method can be used in order to complete the 
final analysis. Different types of MCDA methods are named as 
AHP, ANP, MAUT, MACBETH, PROMETHEE, ELECTRE, 
TOPSIS etc. These are the most commonly used methods for 
multi-criteria decision analysis. The decision-making problems 
can be of different types like, wise ranking problem, rating 
problem, case problem, etc. But in this analysis only rating and 
ranking type analysis is done. So, for this research AHP and 
PROMETHEE methods are selected to carry out sustainability 
analysis for different types of power plants. 

2. Methodology 

2.1 Multi-criteria decision analysis 

MCDA methods are designed in such a way that it helps 
people deciding of complex problems which have many 
alternatives way of solving. When subjective information are 
loaded in the MCDA system it uses it to solve the problem. The 
decision-maker is in control of the system and supplies the 
necessary information to the system. MCDA has a quite large 
range of applications. Mathematical, social science, informatics, 
management, marketing, and economics problems can be solved 
by the MCDA method. Four types of problems are solved widely 
by using the MCDA method. Ranking problems, sorting 
problems, choice problems, description problems. In this paper 
we will use the MCDA to solve ranking and choice problems. 
For solving choice and ranking problems different types of the 
MCDA method has been designed. Some of them are, 
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Table 1 Different ranking and choice-making processes of 
MCDA 


Method Name 

Full meaning 

AHP 

Analytic Hierarchy Process 

ANP 

Analytic Network Process 

MAUT 

The Multi-Attribute Utility Theory 

MACBETH 

Measuring Attractiveness Through a 
Categorical-based Evaluation 

PROMETHEE 

Technique 

Preference Ranking organization 
method for enrichment evaluation 


In our evaluation of the sustainability of different types of 
powerplants, we have used three methods that are AHP, 
PROMETHEE, and TOPSIS. For these methods, the alternatives 
are first selected. There are a total of seven alternatives that are, 
Coal based plant, Nuclear power plant, Oil based power plant, 
Gas power plant, Solar power plant, Hydro power plant, and 
Wind power plant. The main goal is to find the most sustainable 
power plat for electricity generation. This main goal is achieved 
by defining sustainability into some criteria. Then the criteria are 
further sub-divided into some specific sub-criteria. 

After each criterion has been defined each of them are 
needed to be analyzed by the presence of proper data and method 
and thus a final result can be achieved by MCDA method. 

The AHP and PROMETHEE both methods are independent 
of each other. Only one method is enough to get a suitable result, 
but here to methods are used in order to see how much the result 
in individual methods vary from one another. 



Fig. 1 The analysis of most sustainable powerplant 
problem by MCDA method 


2.1.1 Analytic hierarchy process 

AHP is a method that solves a problem by dividing it into 
some criteria and sub-criteria and analysis it. AHP helps by 
giving a simple solution to any complex problem. AHP divides 
a complex problem into simple parts so by solving and analyzing 
these simple parts together, the whole solution can be acquired. 
The problems are needed to be structured and these structured 
problems are needed to compare pairwise for getting a result [3]. 
A hierarchy is designed to solve the problem. The upper element 
is the final decision that is to be taken, the second level is the 


criteria and the criteria is sub divided into some sub-criteria. The 
sub-criteria can also be sub divided if needed. In our case, we 
have divided our goal sustainability analysis of power plants into 
three criteria and eleven sub-criteria. 

After the criteria and sub-criteria are designated then the 
priority of each alternatives is needed to be calculated by a 
priority matrix. A priority scale has been designated by 
American School and this is used in AHP method [4]. 


Table 2 The scale of AHP 


Degree of priority 

Definition 

1 

Equal importance 

3 

Moderate importance 

5 

Strong importance 

7 

Very strong importance 

9 

Extreme importance 


The AHP uses the pairwise comparison between the 
alternatives. The comparison is done manually and then the data 
is collected in a matrix form. The more the number of alternatives 
the more effort is needed to be given. If there are n alternatives 
the number of comparisons for each criterion is (n 2 -n)/2n 2 , total 
comparison that is available to a matrix. Each evaluation with 
alternative with same alternative is regarded as 1. Each 
alternative is been compared with one another and a comparison 
matrix is been built. If we take in consideration the technology 
criteria and efficiency coefficient sub-criteria for AHP analysis 
first the data is collected, here for sample analysis the efficiency 
coefficient sub-criterion has been taken into account. 

2.1.2 PROMETHEE 

The PROMETHEE is also an alternative based choice 
making method that uses data for analysis and is less manual than 
AHP. This method has three basic steps: 

1. Preference degree calculation for each alternative for each 
criterion and sub-criterion; 

2. The uni-criterion flow analysis; 

3. The global flow determination. 

For the first step of calculation the data should be unitarized. 
For this purpose the maximum and minimum value in each data 
set is to be determined and for beneficial criteria, 

x(i)-min{i ) 

X' —— - 

1 max(i)-min(i ) 

For non-beneficial criterion, 

max{i)-x{i ) 

X' —— - 

1 max(i)-min(i ) 

where, 

x t = Value of any alternative for any sub-criterion 

max(Q = Maximum value of an alternative in that sub¬ 
criterion 

min(Q = Minimum value of an alternative in that sub¬ 
criterion 

Then the values each alternatives value in a criterion should 
be subtracted with each other to get a value and then the negative 
value are been neglected. Then each sub-criterion should be 
divided with its own weight. After that each alternative are 
summed and value of the global criterion preference matrix is 
got. 
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Next, the positive and negative flow is calculated. As it is 
hard to draw conclusions from preference degrees the criterion 
preference degree is summarized in a way that is a positive flow 
or noncriterion leaving and negative flow or noncriterion 
entering and the net flows are calculated for each alternative. The 
noncriterion positive flows indicates how an alternative is 
preferred over all other alternatives. The higher the positive flow 
is the more preferred the alternative is. On the other hand, the 
negative flows indicate the average behavior of all alternatives. 
They show how other alternatives are preferred than this 
alternative. By considering both positive and negative flow net 
flow is calculated. The net score lies between -1 to 1. There is 
two ranking method of PROMETHEE in this study 
PROMETHEE-II has been used as PROMETHEE-I doesn’t give 
any comparison value [5]. 

2.1.2.1 The gaia plane 

The Gaia plane is an imaginary plane that contains aspects 
of a decision problem in two dimensional planes. In this Gaia 
plane, actions are represented by bullets and criteria by arrows. 
The position of the actions gives the decision-maker some idea 
about their similarities. The relative position of criterion 
indicates the correlation and anti-correlation between criterions. 
The closer are more correlated. The Gaia plane helps decision¬ 
maker to visualize conflicting point in the process. The length of 
the criterion indicates the power of the generated data. So all the 
information can be got from the Gaia plane in the PROMETHEE 
method [6]. 

3. Data 

All the data are collected from different free sources 
available. Here is the data for different sub-criteria: 

3.1 Technology 

For powerplant technological factors are much important to 
take into consideration. Technology criterion is sub-divided into 
four sub-criteria. They are efficiency co-efficient, energy 
production rate, waste generation and capacity. 


Table 3 Efficiency coefficient of 7 types of power plant [7]. 


Types of powerplant 

Efficiency coefficient 

Coal 

39.4 

Nuclear 

33.5 

Oil 

37.5 

Gas 

43.4 

Hydro 

80.0 

Solar 

9.40 

Wind 

35.0 


3.1.1 Efficiency co-efficient 


The efficiency coefficient can be defined as the ratio of the 
energy output to input energy. Efficiency which has relatively 
high useful energy can be improved by high plant reliability for 
economic benefit. Understandably, high efficiency is consistent 
with high plant reliability and low cost of electricity is 
economically beneficial for powerplant management [8]. As no 
machine cannot give a unit efficiency, there will be losses in 
every powerplant. Here presented a statistical analysis done by 
IAEA 2002 that shows different powerplants efficiency [2]. It 
shows that hydropower has the most efficiency co-efficient 


amongst all the sources and solar has the least. Coal has more 
efficiency co-efficient than nuclear and oil but gas has an upper 
hand on coal. 

3.1.2 Energy production 

Energy production rate can be defined as the amount of 
energy produced by per kg of fuel burning in per hour time. The 
calorimetric analysis is used to determine the energy production 
rate of any fuel [9]. A report published by the University of 
Washington, gives data on energy production of different sources 
if other factors are constant. It can be seen that gas has the more 
energy production rate than the other sources and can be more 
beneficial in selection. Nuclear is more efficient than the other 
non-renewable and the renewable sources have infinity amount 
of energy production. 

Table 4 Energy production of different power plant [9]. 


Types of power plant Energy production (KW.hr/lb) 


Coal 

4.4 

Nuclear 

10.6*10 A 3 

Oil 

6.67 

Gas 

19.07 

Hydro 

Infinity 

Solar 

Infinity 

Wind 

Infinity 

3.1.3 Waste generation 


In this evaluation, the waste generation is taken as the C02 
emission from the different power plants. Measuring the life- 
cycle greenhouse gas emissions involves measuring the global 
warming potential of electrical energy sources by assessing the 
life-cycle of each source. The findings are described in units of 
potential for global warming per unit of produced electrical 
energy from that source. Waste generation can be expressed in 

gram equivalent of C02 in per 

KW of energy production. 

Thomas Bruckner in an article called “Mitigation of Climate 

Change” published the for C02 

emission from the different 

power plants [10]. Coal powerplants emit the most amount of 
carbon-die-oxide and nuclear and wind are the cleanest source of 

energy. 


Table 5 Waste generation from different powerplants [10]. 

Types of power plant 

Waste generation 
(gC02eq / kWh) 

Coal 

820 

Nuclear 

12 

Oil 

520 

Gas 

490 

Hydro 

24 

Solar 

48 

Wind 

12 


3.1.4 Capacity factor 


The capacity factor is a dimensionless parameter used in the 
evaluation of power plants. The capacity factor of a powerplant 
can be defined as the actual amount of energy produced in a 
certain interval of time and the maximum energy that could have 
been produced by running the power plant at full power at that 
time. This parameter or criteria has great importance in power 
plant evaluation. The capacity generally is less than 100%, as 
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machines need to do periodic maintenance and sometimes 
corrective maintenance too. The performance of a power plant 
greatly depends on its capacity. It is a beneficial criterion. On a 
report published by El A, 2019 different powerplant capacity was 
analyzed [11]. It has been found that nuclear power plant has the 
height capacity and the solar has the lowest. 


Table 6 Capacity factor of different power plant [11]. 


Types of power plant 

Capacity factor 

Coal 

85 

Nuclear 

90 

Oil 

30 

Gas 

87 

Hydro 

50 

Solar 

20 

Wind 

38 


3.2 Safety and sustainability 


In the evaluation of powerplants, the safety of human life 
and the plant equipment must be taken into consideration. Also, 
the sustainability of energy sources should be ensured. The sub¬ 
criteria that are been in this criterion is availability, reserve to 
production ratio and accidental death analysis. 

3.2.1 Availability 

In the case of powerplant, the availability describes the time 
period that a power plant is able to generate electricity over a 
certain period of time. A power plant can be out of service due 
to maintenance or repair [12] or due to weather conditions like 
lack of sunlight or wind etc. [13]. The equipment quality, 
maintenance, types of fuel used and how the operation is done in 
the power plant greatly affects the availability of the power plant. 
If other factors are equal power plants which run more do have a 
higher value of availability. IAEA published a report in 2002 
[14] that shows availability of various types power plants. Here 
it can be seen that nuclear has the most availability on the other 
hand solar has the least once. Oil and gas dominate over coal 
power plant. 

Table 7 Availability of different power plant [14]. 


Types of power plant Availability 


Coal 

85.4 

Nuclear 

96 

Oil 

92 

Gas 

91 

Hydro 

50 

Solar 

20 

Wind 

38 


3.2.2 Reserve to production ratio 


The reserve to production ratio is the ratio of the currently 
available amount of fuel to the annual consumption of each non¬ 
renewable resources. When considering the fuel sources only 
well- known sources are taken into account. Different types of 
models are used in this case. Non-economic extraction of fuel is 
not considered as available sources. British petroleum publishes 
a report every year based on reserve to production ratio of every 
available non-renewable sources. Here from the BP report 2019 
[15] all the data has been taken for 7 types of powerplants. Here 
it can be seen that all renewable sources have the maximum 


infinity reserve to production, oil has the most minimum 
resources compare to other. These values represent the current 
situation but new sources of energy are been discovered always. 


Table 8 Rpr for different power plant [15]. 


Types of power plant 

RPR 

Coal 

164 

Nuclear 

70 

Oil 

40.5 

Gas 

66.7 

Hydro 

Infinity 

Solar 

Infinity 

Wind 

Infinity 

3.2.3 Accidental analysis 


One of the vital concerns in the power plant is to ensure the 
safety of the workers and the public. For this, the powerplant with 
the least amount of accidental deaths is preferable as human life 
is more superior to power generation. For the accidental analysis 
purpose, total amounts of deaths are taken into account until now 
from different power sources in different types of accidents. 
Organization for economic and co-operation development 
published a report that shows deaths from different power plants 
[16] The maximum amount of deaths are caused by coal-based 

power plant. 


Table 9 Accidental deaths from different power plants [13]. 

Types of power plant 

Accidental Death 

Coal 

36441 

Nuclear 

31 

Oil 

20218 

Gas 

2043 

Hydro 

29938 

Solar 

0 

Wind 

0 

3.3 Economy 


For business purposes, the 

economy should be taken into 

concern for powerplants as the 

main purpose of a powerplant 

should produce energy at a lower cost. So for selecting the best 
powerplant all costs should be taken into account. The costs that 
are important in the powerplant economy are capital costs, the 
operation, and management cost, the fuel cost, the external cost, 

etc. 


Table 10 Capital costs of different power plants [17]. 

Types of power plant 

Capital Costs (€/KW) 

Coal 

3661 

Nuclear 

6016 

Oil 

1802 

Gas 

1079 

Hydro 

2752 

Solar 

7191 

Wind 

5446 


3.3.1 Capital cost 


The capital cost can be regarded as the summation of land 
cost, the necessary cost for infrastructure and the necessary 
equipment costs that are needed to run the power plant. Labor 
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costs needed to build the powerplant are not included here. EIA 
published a report in 2019, publishing the total capital cost 
needed to build the powerplant [17]. It can be seen that coal and 
nuclear have relatively high capital costs than others. Oil power 
plants have lowest capital costs among other electricity 
generation options. 

3.3.2 Operations and maintenance cost 

Wages of employees, energy funds, the products and service 
for power plant operation include the O&M cost. Also, the 
management cost is included in this cost as it is very much 
necessary for proper maintenance for the prolonged 
serviceability of the power plant. Also, if any part of the power 
plant gets deuterated it is needed to be fixed. The operation 
management cost is divided into two sub categories; fixed O&M 
cost and variable O&M cost. Fixed operation cost is the cost that 
is calculated yearly based and has no effect on the amount of 
energy produced and the variable cost is directly related to the 
amount of energy produced. The value of O&M cost is got from 
the IAE report 2019 [17]. The variable O&M cost are zero for 
renewable sources as while production they don’t need extra 
maintenance. 


Table 11 Fixed O&M cost of different power plants [17]. 


Types of power plant 

O&M Cost (€/KWyr) 

Coal 

4.613 

Nuclear 

13.83 

Oil 

4 

Gas 

1.6 

Hydro 

4.75 

Solar 

9.7 

Wind 

12.5 

Table 12 Variable O&M cost of power plants [17]. 

Types of power plant 

O&M Cost (€/KWyr) 

Coal 

4.48 

Nuclear 

2.36 

Oil 

5.67 

Gas 

2.54 

Hydro 

1.39 

Solar 

0 

Wind 

0 

Table 13 Fuel cost of different power plants [14]. 

Types of power plant 

Fuel Cost (€cent/KWyr) 

Coal 

1.31 

Nuclear 

0.27 

Oil 

1.84 

Gas 

2.34 

Hydro 

0 

Solar 

0 

Wind 

0 


3.3.3 Fuel cost 


Fuel is something that is burned and produces energy. For 
power plant economy it is always required the cheapest value 
fuel. The cheaper the fuel is the more economic benefit the plant 
has. Fuel releases energy from the chemical reaction and passes 
this energy to run the turbine and produces electric energy. The 
discharge is performed in such a controlled way that this energy 


can utilize to achieve maximum work. IAEA report 2002 gives a 
basic price value of all fuels [14]. It can be seen that renewable 
sources have less fuel cost and gas is the most expensive one. 

3.3.4 External cost 

These costs are categorized with the extra cost of health and 
the environment. It is not mixed with the cost of the generation 
of electricity. These external costs are related to the restoration 
of the negative side effects that happen on human health as well 
as the ecosystem during the operation power plant. They are 
calculated on the basis of the life cycle external cost of the 
powerplant. By report IAEA 2002 an overview data of different 
powerplants is got [14]. 


Table 14 External cost of different power plants [14]. 


Types of power plant 

External Cost(€cent/KWyr) 

Coal 

8.4 

Nuclear 

0.49 

Oil 

6.75 

Gas 

2 

Hydro 

0.56 

Solar 

0.24 

Wind 

0.16 


4. Results and discussion 


Here in this analysis two cases have been considered. Case 
1 Global criterion and sub-criterion mean taking each criterion 
with equal importance and the relatively analyzing the sub¬ 
criterion based on their relative importance. And case 2 means 
taking all sub criterion as equal. These cases are not constant, 
they are both independent. This case varies with time, location 
etc. So as for sample analysis two basic cases have been 
considered. Either one of the case can give suitable result. 

4.1 Analysis conditions and results by AHP 

For the analysis purpose different sample condition were 
taken by varying the weight of criterion. As by multiplying the 
weight factor with the priority factor the importance of each 
alternatives can be got. For analysis purpose two conditions are 
taken for every method. 

4.1.1 Case 1 (Global sub-criterion) 

For this condition technology, safety & sustainability and 
economy all the criteria have been considered with an equal 
weight of 33.33%. Next for determining the weight of each sub¬ 
criterion pairwise comparison is done. If the economy criterion 
is taken as an example every sub-criterion is compared with other 
pairwise and the pairwise value is got. 


Table 15 Pairwise comparison of economy criteria 


Sub criterion 
of Economy 

Capital 

cost 

Operation & 
management 
costs 

Fuel 

costs 

External 

costs 

Capital cost 

1 

5 

1 

3 

Operation & 
management 

1/5 

1 

1/5 

1/3 

costs 

Fuel costs 

1 

5 

1 

3 

External 

Costs 

1/3 

3 

1/3 

1 
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This is how every criterion is analyzed and a global sub- 
criteria index is got by the process. As the pairwise comparison 
is easy and no exact data are not needed for determining the 
weight, this method is vastly used in the process of determining 
the weights. The global criteria mainly give information about 
these three equal criteria. A global ranking of 100% of power 
plant evaluation different weight of sub-criteria is got which is 
shown on Table 16. 

Table 16 Global Criterion and sub-criterion weight 


Sustainability of power plant (100%) 


Technology 

(33.33%) 

Safety and 
sustainability 
(33.33%) 

Economy (33.33%) 

Energy 

coefficient 

2.33% 

Availability 

4.30% 

Capital 

Costs 

12.92% 

Waste 

Generation 

17.5% 

Reserve to 
production 
ratio 

24.71% 

Operation & 
Management 
Costs 

2.20% 

Energy 

Production 

9.00% 

Accidental 

death 

4.30% 

Fuel Costs 

12.92% 

Capacity 

4.55% 

analysis 


External 

Costs 

5.91% 


By applying global criterion and sub-criterion weight for 
different powerplants the final results are got by AHP analysis. 
It can be seen that the wind power plant is more preferable than 
the other power plants and hydro is also same preferability like 
wind. Oil is less preferable than the other sources and nuclear has 
an advantage on coal. 


11 | Overall Percentage score | 

20.51 



Coal Nuclear Oil Gas Hydro Solar Wind 
Power plant 


Fig. 2 AHP analysis by Case 1 Global sub-criterion (Equal 
criteria) 

4.1.2 Case 2 (All sub-criteria is of the same weight) 

In this case of analysis, all the sub-criteria’s weight is 
regarded as same and they are 9.09% each for every eleven sub¬ 
criteria. The result shows that hydro is the most beneficial source 
if all the sub-criteria are taken as equal and gas has more 
advantage than other renewable sources like oil, coal 


I | Overall Percentage score | 




Coal Nuclear Oil Gas Hydro Solar Wind 
Power plant 


Fig. 3 AHP analysis taking each sub-criterion equal 
4.2 Analysis conditions and results by PROMETHEE 

In PROMETHEE method the result is shown on the Gaia 
plane with all scatter points. The higher the position in the y axis 
the more preferable the power plant is. The Gaia plane is a two- 
dimensional plane that shows different point in x-y axis with 
preferable value. 

4.2.1 Case 1 (Global sub-criteria weight for every criterion has 
equal importance) 

The value of global flow was determined for global sub¬ 
criteria weight by taking each criterion as 33.33%. By the 
PROMETHEE method analysis, it can be seen that hydropower 
has the most global flow than other sources. Renewable sources 
have always the great acceptancy than other sources. Oil has the 
least amount of acceptancy and nuclear has more acceptancy 
than coal and oil. 
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Power plant 


Fig. 3 PROMETHEE analysis by global sub-criteria (Equal 
criteria) 

4.2.2 Case 2 (All sub-criteria is of same weight) 

As for this every sub-criterion is 9.09% the global flows 
show the priority of global flow index for 7 types of power plant. 
It can be seen that hydro is more efficient than other sources. 
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Fig. 4 PROMETHEE analysis taking each sub-criterion 
equal 

4.3 Basic comparison between the AHP and PROMETHEE 
result 

In both AHP and PROMETHEE analysis it has been seen 
that the renewable sources are better than non-renewable sources 
and also nuclear has far better sustainability than other sources 
of non-renewable. From the average result analysis, it can be 
seen that hydro is the more sustainable source to build and oil is 
less sustainable source to build. Nuclear power plant plants are 
much better in case of sustainability. 

5. Conclusion 

Here the sustainability part was defined properly in the case 
of power plant and the sustainability of a power plant is divided 
into some criteria and those-criteria is further sub-divided into 
several sub-criteria and by this, a hierarchy tree was made. 
Analytical Hierarchy method and PROMETHEE was studied 
properly and finally the sustainability index was made. By the 
analysis of seven types of power plants with respect to 
technology, economy, and safety & sustainability with the help 
of the multi-criteria decision analysis for several cases, it can be 
seen that renewable sources are more sustainable at present. 
Renewable sources have an endless supply and they can serve as 
a spontaneous source of energy. Hydro, solar and wind are the 
highest-ranked power plant in terms of sustainability. The 
nuclear power plant is also a very beneficial in terms of 
sustainability and some cases is even more beneficial than 
renewable sources. Coal has the greatest reserve amongst all the 
non-renewable fuels analyzed so it can be advantageous too. Gas 
has an upper hand on coal due to its low capital cost. 
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ABSTRACT 


The performance of diesel engine highly depends on atomization, vaporization and mixing of fuel with air. These factors are strongly 
influenced by various parameters e.g. injection pressure, injection timing, compression ratio, equivalence ratio, cylinder geometry, in 
cylinder air motion etc. In this study, a diesel engine has been investigated by employing a commercial CFD software (ANSYS Forte, 
version 18.1) especially developed for internal combustion engines (ICE) modeling; focusing primarily on the effects of equivalence 
ratio and compression ratio on combustion and emission characteristics. RNG k-s model was employed as the turbulence model for 
analyzing the physical phenomena involved in the change of kinetic energy. In order to reduce the computational cost and time, a sector 
mesh of 45o angle with periodic boundary conditions applied at the periodic faces of the sector, is considered instead of using the whole 
engine geometry. Simulations are performed for a range of equivalence ratio varying from 0.6 to 1.2 and for three compression ratios 
namely, 15:1, 18:1 and 21:1. Results show that, improvement in combustion characteristics with higher compression ratio could be 
achieved for both lean and rich mixtures. Peak in-cylinder pressure and peak heat release nearer to TDC are achieved for compression 
ratio of 18:1 that could results in more engine torque. For compression ratio beyond 16:1, effects of fuel concentration on ignition delay 
is more pronounced. At lower compression ratio, in-cylinder temperature is not sufficiently high for atomization, vaporization, mixing 
of fuel with air, and preflame reactions to occur immediately after the fuel injection. NOx emission in diesel engine increases due to 
higher pressure and temperature inside the cylinder associated with relatively higher compression ratio. Rich mixture leads to more CO 
and unburnt hydrocarbon emission compared to lean mixture as result of incomplete combustion. Engine operation with too high 
compression ratio is detrimental as emission is a major concern. 


Keywords: Computational Fluid Dynamics (CFD), Diesel engine, Compression ratio, Combustion performance, Engine 
emissions. 
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1. Introduction 

Due to limited reserve of petroleum, ever increasing energy 
demand and stringent emissions legislation; automakers are 
striving to develop advanced engine with improved the fuel 
economy, reduced emissions and higher thermal efficiency. 
Diesel engine is a high compression ratio internal combustion 
engine where auto-ignition of fuel occurs due to high 
temperature developed in the compression stroke. Diesel engine 
combustion is considered as the series of complex processes and 
the main events are categorized into ignition delay, uncontrolled 
rapid combustion phase, mixing controlled combustion phase 
and after burning phase. The time interval between start of fuel 
injection and start of combustion refers to ignition delay. This 
delay period can further be divided into two parts - Physical 
delay and chemical delay. The fuel, accumulated during physical 
and chemical delay, is burnt rapidly causing very high heat 
release in uncontrolled combustion phase. After that the burning 
rate is controlled by the rate of fuel injection and mixing of fuel 
with the surrounding air, in controlled combustion phase. After 
burning phase involves lower rate of heat release due to 
combustion of remaining fuel [1]. In transportation sector, diesel 
engine is preferred over petrol engine because of its higher 
thermal efficiency, better fuel economy and higher engine power 
and torque. However, they suffer from excessive emissions 
which require the use of expensive exhaust-gas after-treatment 
devices [2]. The emission products include soot particles, oxides 


of nitrogen (NOx), carbon monoxides (CO), unbumt 
hydrocarbons (UHC) and particulate matters (PM) [3]-[6]. Due 
to the incomplete combustion of fuel unburnt hydrocarbon 
(UHC) and carbon monoxides (CO) are generated. The amount 
of NOx particles, soot products and other undesired gases in the 
exhaust gas could be higher than the gasoline engine. 

In recent years, extensive researches on diesel engines are 
ongoing aiming to reduce NOx and soot emission without 
compromising higher thermal efficiency of the engine. In a study 
by Reitz et al. [7], the effects of EGR, equivalence ratio and 
mixing timing on emission levels for heavy duty Premixed 
Charge Compression Ignition (PCCI) engine has been analysed. 
At high EGR percentage, drastic reduction in PM levels is 
reported. Intake temperature is increased with high EGR rate 
which reduces the PM level. The emissions of PM, NOx and HC 
are reduced and emission of CO is increased with the increase of 
mixing timing. As mixing time increases, the charge (fuel-air) 
becomes more homogeneous and thus reduces emissions. The 
emission of PM is the function of in-cylinder average 
equivalence ratio. It is seen that at lean equivalence ratio, PM is 
not formed. Jindal et al. [8] investigated the effect of 
compression ratio and injection pressure in a direct injection 
diesel engine running on jatropha methyl ester. Being a fuel of 
different origin, the standard design parameters of diesel engine 
may not be suitable. For small size direct injection constant 
pressure engine, optimum combination of emission is found at 
CR 18 with injection pressure of 250 bars. Raheman et al. [9] 
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studied the performance of Richardo E6 engine using biodiesel 
obtained from mahua oil (B100) and its blend with high speed 
diesel (HSD) at varying compression ratio (CR), injection timing 
and engine loading. Brake specific fuel consumption (bsfc) and 
exhaust gas temperature increases, whereas brake thermal 
efficiency decreases with the increase in proportion of biodiesel 
in the blends at all compression ratios (18:1 to 20:1) and injection 
timings (35-45° before TDC). Effect of equivalence ratio on 
combustion and emissions with a low compression ratio of 14.2 
in a dual-fual 6-cylinder engine was investigated experimentally 
by Jinbao et al. [10]. They reported that the gas consumption rate 
decreases with the increase of equivalence ratio regardless of 
strategy in nozzle parameter, exhaust gas recirculation (EGR) 
rate, injection parameter. Peak heat release rate and exhaust gas 
temperature increased whereas combustion duration decreased 
due to the high equivalence ratio in dual fuel ignition mode. Juan 
et al. [11] studied the knocking tendency in a Cooperative Fuel 
Research (CFR) engine due to the effect of the varying 
equivalence ratio. They reported that lean mixture reduces the 
knocking tendency allowing to operate the CFR engine at higher 
critical compression ratio. Dong Nam [12] carried out post¬ 
graduate research on reducing engine emissions (NOx) from 
marine diesel engines. He concluded that NOx emission can be 
reduced by primary methods such as retard injection, fuel nozzle 
modification, change of compression ratio, water direct injection, 
water emulsification, exhaust gas recirculation (EGR) and 
secondary method such as selective catalytic reduction (SCR). 
Jafarmadar [13] worked on NOx reduction by split injection. In 
diesel engine, to gain reduction in NOx emission without 
increasing soot emission is very difficult by adjusting injection 
timing only. This study demonstrated that split injection could be 
a promising technique to reduce both NOx and soot emission 
simultaneously. Recently, numerical simulation of ICE through 
Computational Fluid Dynamics (CFD) modelling, has become a 
promising tool to get better insight about the complex processes 
such as fuel injection, spray formation, chemical kinetics, 
ignition model, flame propagation, knock phenomena, emission 
formations that occur inside the cylinder of ICE. Taqizadeh et al. 
[14] performed numerical investigation of effects of equivalence 
ratio on combustion characteristics of an RCCI engine. They 
found that, by increasing equivalence ratio from 0.35 to 0.55 in 
a constant energy ratio, noticeable growth in the maximum 
amount of pressure and temperature could be achieved; 
consequently, NOx emission would increase significantly, IMEP 
increases by 43%, and ISFC decreases by 30%. 

Though many experimental investigations are reported so 
far, but only a few computational fluid dynamic studies on diesel 
engines are available in the literature. Therefore, in the present 
study, effects of equivalence ratio and compression ratio on 
combustion and emission characteristics in a diesel engine have 
been investigated through computational fluid dynamics (CFD) 
modelling. 

2. Computational approach 

ANSYS FORTE, version 18.1, has been used in this 
numerical study of diesel engine operating with varying 
equivalence ratio and compression ratio. This is a fluid-dynamics 
software that specializes in the simulation of combustion 
processes in an internal combustion engine, using a highly 
efficient coupling of detailed chemical kinetics, liquid fuel spray 
and turbulent gas dynamics. The software solves the full 
Reynolds-averaged Navier-Stokes (RANS) equations with well- 


established flow turbulence models. RNG k-s model was 
employed as the turbulence model for analyzing the physical 
phenomena involved in the change of kinetic energy. 

2.1 Basic governing equations 

The working fluids for internal combustion engines are air- 
fuel mixture before combustion and the burned products after 
combustion. Basic fluid dynamics are governed by the Navier- 
Stokes equation in the turbulent reacting flows. Governing 
equations are mainly derived by- 

• Assuming ideal gas law for gas phase 

• The use of Fick’s law for mass diffusion 

• Assuming Newtonian fluid 

• The use of Fourier’s law for thermal diffusion 

2.1.1 Species conservation equation 

dp k 

— + V(p k u) = V[p k Dy k ] + Vcp +p c k 

+ p s k (Here,k = 1,2, ...K) 

Here, p is the density, subscript k means species index, K is 
the total number of species, u is the flow velocity vector 

The term cp accounts for the effects for the ensemble¬ 
averaging the convection term i.e. 


cp = p k u- p k u. 

p k and p k are terms due to chemical reactions and spray 
evaporation, respectively. 

2.1.2 Fluid continuity equation 

The continuity equation fir the total gas phase fluid is, 

dp _ 

— + V(pu) = p s 

2.1.3 Momentum conservation equation 


The momentum equation considering the effects of 
convection, pressure force, viscous stress and turbulent transport 
as well as the impact from liquid sprays and body force is, 


dpu _ 2 _ __ 

——h V(puu) = -Vp + Vo — -pkl + F s 4- pg 
ot 3 

where, p is the pressure, F s is the the rate of momentum gain 
per unit volume due to the spray, g is the specific body force, a 

is the viscous shear stress, <7 = pv [vfi + (Vfi) r — ^ (Vfi)/j 

V= laminar kinetic velocity, 1= identity tensor, superscript T 
means transpose of a tensor. 

The stress accounts for the effects for ensemble averaging 
of filtering of the nonlinear convection term i.e. F 1 = p(wu — 
uu). In the RANS approach it is called Reynolds Stress. 


2.1.4 Energy conservation equation 

Based on First Law of Thermodynamics, the change of 
internal energy has to be balanced by the pressure work and heat 
transfer. The internal energy transport equation is. 

dpi , __ ~ - r 

——h V(p ul ) = —pV u — pVu — V] - VH 4- pi 4- Q c 
ot 

+ Q S 
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where, I is the specific internal energy, J is the heat flux 
vector, J = —AVf - pD T h k V ( j) 

A is the thermal conductivity, which is related to the thermal 
diffusivity a and heat capacity C p by A = pC p a, T is the 
temperature and h k is the specific enthalpy of species k. £ is the 
dissipation rate due to chemical heat release and spray 
interactions, respectively. 

The term H accounts for the effects of ensemble averaging 
or filtering of the convection term, i.e. H = p(ul — ul ) 

2.1.5 Reynolds averaged navier stokes (RANS) approach 


The aim of the RANS is to simulate ensemble-average flow 
field. The ensemble average is an average that allows the 
interpretation of the phenomena in terms of the repeatability of 
multi-component flows. 

To model the turbulent transport process gradient-diffusion 
assumptions are approached. For momentum equation, the 
deviatoric components of the Reynolds stress are assumed to be 
proportional to the mean deviatoric rate of strain. The Reynolds 
stress tensor is, 


f = —pd T 


Vu(Vu) T — -{Vu)I 


2 

+-P ki 


where, d T is the turbulent kinematic viscosity, if dissipation 
rate is £ then, 


P 



k is the turbulent kinetic energy, defined by, 


k = — trace (T) = -uu 

2 p 2 

C p is a model constant that varies in different turbulence 
model formulations. 

The turbulent flux term in the species transport equation is 
modeled as, (p = pD T Vy k , where D T is the turbulent diffusivity. 

Similarly, turbulent flux term 

H = —AVf - pD T ^ h k Vy k 

k 

where y k is the turbulent conductivity and is related to the 
turbulent thermal diffusivity a T and heat capacity C p by A T = 
pC p a T 

The turbulent mass and thermal diffusivity are related to the 
turbulent viscosity by, 


D r = 


■>c T 


d T 

a T = y- 
r r T 


where, 5^ and P rT are the turbulent Schmidt and Prandtl 
numbers, respectively. 

2.2 Geometry and mesh generation 


In case of diesel engine with equally spaced multi-hole injector, 
properties such as injection pressure, temperature, air-fuel 


mixture homogeneity or mixture concentration, combustion 
process etc. are considered identical with respect to each hole in 
the injector and the corresponding spray. 



Fig. 1 A 45° sector mesh for modeling diesel engine. 

Therefore, the combustion chamber inside the engine 
cylinder is thought to be symmetrical with respect to each nozzle- 
hole in the injector. As a consequence, an advantage can be 
readily obtained by simulating a certain section of geometry 
instead of considering whole engine geometry for calculation. 
For example, an eight-hole injector allows simulation using a 45° 
sector (36078). In the current study, for computational domain 
and mesh generation, a sector of 45° angle with periodic 
boundary conditions applied at the periodic faces of the sector is 
considered instead of using the whole geometry of the engine in 
order to reduce computational cost without comprising the 
accuracy (Fig. 1). 


Table 1 Engine specifications and operating conditions 


Engine base type 

Cummins N-14, DI diesel 

Number of cylinders 

1 

Cycle 

4-stroke 

Number of intake valves 

2 

Number of exhaust valves 

1 

Combustion chamber 

Quiescent, direct injection 

Squish 

0.56 cm 

Bore 

13.97cm 

Stroke 

16.02 (CR=15), 16.21 
(CR=18), 16.35 (CR=21) 

Range of compression 
ratio 

15:1, 18:1,21:1 

Bowl width x depth 

9.79 x 1.55 cm 

Connecting rod length 

30.2 cm 

Engine Speed 

1200 rpm 

Range of equivalence ratio 

0.6- 1.2 

Fuel injector type 

Common rail, pilot valve 
actuated 

Number of holes 

8, equally spaced 

Nozzle orifice area 

3.02x1 O' 4 cm 2 

Start of Injection (SOI) 

-22° [After TDC] 

Inflow Droplet 

Temperature 

368.0 K 

Discharge Coefficient 

0.7 
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A single cylinder DI diesel engine based on a Cummins N- 
series production engine has been modelled by using ANSYS 
Forte CFD simulation software. Details of the engine is available 
in Refs. [15],[16]. The engine specifications and operating 
conditions are given in Table 1 . 

2.3 Fuel chemistry and sub-models 

For specifying the chemical reaction mechanism for fuel 
combustion ANSYS Forte supports the industry standard 
CHEMKIN format. A reduced n-heptane mechanism that can be 
used to represent the diesel fuel under conventional diesel-engine 
combustion conditions and available in fuel library of ANSYS 
Forte is added to the chemistry model in the CFD code to 
simulate the combustion of diesel fuel. Some other important 
CFD sub models used in this simulation work are shown in Table 
2 . 


Table 2 CFD Sub-models 


Turbulence model 

RNG k — a model 

Droplet breakup model 

KH-RT model 

Combustion model 

Mixing controlled combustion 
(MCC) 

Fuel chemistry model 

n-heptane reduced mechanism 

Soot model 

Two-step semi-empirical model 

NOx formation 
mechanism 

Zeldovich mechanism 


2.4 Initial and boundary conditions 


The initial conditions in the engine simulations are presented 
in Table 3. These values are utilized for cases with different 
diesel injection timings and gasoline/diesel ratio. Law of the 
Wall model is specified for piston, head and liner which can 
capture wall boundary layer effects more accurately for 
boundary layers that are thinner than the mesh size. For in¬ 
cylinder engine simulations, the turbulent law-of-the-wall 
velocity condition and fixed temperature walls are usually 
employed. 

Table 3 Initial conditions 


Parameters 

Value 

Intake Valve Closing 

165° before TDC 

Exhaust Valve Open 

125° after TDC 

Temperature at Inlet Valve Closing 

362 K 

Pressure at Inlet Valve Closing 

2.215 bar 

Turbulent Kinetic Energy 

10,000 cm 2 /sec 2 

Turbulent Length Scale (cm) 

1.0 

Initial Swirl Ratio 

0.5 


Table 4 Boundary conditions 

Periodicity 

45 degree 

Wall model 

Law of the wall 

Piston Temperature 

500 K 

Head Temperature 

470 K 

Line Temperature 

420 K 


3. Results and discussion 

3.1 Effect of compression ratio on combustion characteristics: 

For diesel engines, the engine combustion is often simulated 
from intake valve closure (IVC) to exhaust valve opening 
(EVO), rather than modeling the full air intake or exhaust flow 
processes involving the intake and exhaust ports, respectively. 
This is usually a reasonable approximation since the gas in the 
cylinder at IVC is a relatively homogeneous mixture of air and 
exhaust gas (due to internal residual or from exhaust-gas 
recycling), prior to fuel injection. Here, in this study only a part 
of the engine operation, namely, in between IVC (i.e. 165 degree 
before TDC) to EVO (i.e. 125 degree after TDC) is simulated for 
the desired investigations. 

Fig. 2 show the combustion characteristics of diesel in terms 
of in-cylinder pressure, temperature and apparent heat release 
rate (AHH) for compression ratio of 15, 18 and 21. Both lean 
mixture with equivalence ratio, (p =0.8 and rich mixture with 
equivalence ratio, (p =1.2 are considered to analyse the effect of 
fuel concentration. It is observed that peak in-cylinder pressure, 
pressure rise rate, peak in-cylinder temperature, rate of 
temperature rise and heat release rate, all increase with 
compression ratio for both lean and rich mixtures. It is well 
known that, in Cl engine prior to fuel injection both air 
temperature inside the cylinder and cylinder wall (including 
head, liner and piston top) temperature will be greater for higher 
compression ratio. This high temperature increases the mixture 
reactivity, and rate of mixing of fuel with air, causing rapid 
combustion of fuel that lead to steeper heat release rate, in¬ 
cylinder pressure and temperature. Among the three compression 
ratios, CR=15, CR=18 and CR=21; both peak in-cylinder 
pressure and peak in heat release are close to TDC for CR=18. It 
is desirable to have both peak in-cylinder pressure and peak in 
heat release nearer or closes to TDC as much as possible because 
these characteristics will lead to maximum engine power or 
torque. Therefore, it could be concluded that, in terms of 
combustion characteristics, compression ratio, CR=18 is the 
optimum one. 


Initial Swirl Profile Factor 3.11 
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Fig. 2 In-cylinder Pressure, Temperature and Apparent Heat Release Rate (AHRR) Vs Crank angle at different equivalence ratio 

and compression ratio (for CR=15, CR=18 and CR=21) 
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After start of injection (SOI) of fuel, the time delay up to 
start the combustion (SOC) is termed as ignition delay which can 
be subdivided into physical delay and chemical delay. The higher 
temperature due to high compression ratio will facilitate 
atomization, vaporization and mixing of fuel with entraining air 
inside the combustion chamber, thereby decreasing physical 
delay. On the other hand, chemical delay also become shorter as 
high temperature accelerate the pre-flame reactions and make the 
air-fuel mixture ready for autoignition. Ignition delay becomes 
shorter with compression ratio over a range of equivalence ratio 
as depicted in Fig. 3. For compression ratio beyond CR=16, rich 
mixture shows lower ignition delay compared to lean mixture. 
At lower compression ratio, for example CR=15, it takes around 
17-19 deg. of crank angle (CA) after SOI at 22.5° before top 
dead center, TDC (0° crank angle represents TDC). This is 
because in-cylinder temperature is not sufficiently high for 
atomization, vaporization, mixing of fuel with air, and preflame 
reactions to occur immediately after the fuel injection. 



Fig. 3 Ignition delay as function of compression ratio for 
different equivalence ratios 

3.2 Effect of compression ratio on emission characteristics: 

It is evident from Fig. 4 that diesel engine NOx emission 
increases for relatively higher compression ratio. Previously it is 
seen that higher compression ratio causes higher pressure and 
temperature inside the cylinder. At high temperature, the 
nitrogen of the air can react with oxygen and forms various 


nitrogen oxides e.g. nitric oxide (NO) and nitrogen dioxide 
(N02). According to thermal mechanism (commonly known as 
Zeldovich mechanism) temperature and oxygen concentration 
are the dominant factors in NOx formation; and it is found to 
increase exponentially with in-cylinder temperature [1]. This 
mechanism is highly sensitive to temperature not only because 
of the high activation energy of reaction but also as a result of the 
rapidly increasing oxygen atoms concentration in flames with 
increasing temperature [17]. Fig. 4 clearly demonstrates that for 
lean mixture with equivalence ratio, (p =0.8, the exhaust emission 
of NOx is relatively higher but combustion of rich mixture ((p 
=1.2), reduces NOx emission drastically. As combustion is 
started earlier with higher compression ratio, therefore, emission 
of CO is seen to emit earlier in case of CR= 21 and CR=18 
compared to CR=15 as shown in the Fig. 4. Besides, introducing 
relatively higher amount of diesel fuel (0 =1.2), results in more 
CO emission compared to lean mixture with equivalence ratio, (p 
=0.8 which indicates incomplete combustion is more pronounced 
in case of combustion of rich mixture. Unburnt hydrocarbon 
emission is observed more in lower compression ratio, CR=15 as 
a result of poor oxidation associated with lower in-cylinder 
temperature. 

3.3 Effect of equivalence ratio on combustion characteristics: 

For diesel engines, fuel concentration plays a vital role on 
combustion performance in terms of in-cylinder pressure, 
temperature and apparent heat release rate (AHRR). Fig. 5 
illustrates that for fuel rich mixture, combustion of diesel fuel 
results in faster rate of pressure rise. As more fuel is injected for 
(p =1.0 and 1.2, compared to lean mixture, (p =0.6 to 0.8, 
therefore, a combustible mixture is readily formed in the 
combustion chamber and autoignition of fuel-air mixture occurs 
earlier. Though energy supplied is greater for comparatively rich 
mixture but combustion temperature is found to be lower for (p 
=1.0 and 1.2, compared to lean mixture, (p =0.6 to 0.8. This might 
be attributed to the fact that even after autoignition of fuel-air, 
the injection of diesel fuel continues to occur directly in burning 
zone; as a consequence, a portion of the heat of combustion is 
absorbed by newly injected fuel droplets before vaporization. 
The highest combustion temperature is observed for (p = 0.7 and 
(p = 0.8. In addition, both peak in-cylinder pressure and peak in 
heat release are close to TDC for (p = 0.7 and (p = 0.8, which will 
lead to maximum power or torque. 
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Fig. 4 Emission of NO x , CO and Unburnt Hydrocarbon Vs Crank angle at different equivalence ratio and compression ratio (for 

CR=15, CR=18 and CR=21) 
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Fig. 5 In-cylinder Pressure, Temperature and Apparent Heat Release Rate (AHRR) Vs Crank angle over a range of equivalence 

ratio and compression ratio 
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3.4 Effect of equivalence ratio on emission characteristics: 

As mentioned before NOx formation inside the combustion 
chamber is highly temperature dependent process. It is also well 
known that lean burn strategy in diesel engine leads to high NOx 
emission. In the previous section, the in-cylinder temperature is 
found to be greater for equivalence ratio in the range of 0 = 0.6 
to 0 =0.8. In this numerical modelling, highest level of NOx 
emission is observed for lean combustion (0 = 0.6 and (j) = 0.7), 
as can be seen in Fig. 6. Whereas, for 0 = 1.0 and 0 = 1.2, due to 
low combustion temperature, NOx formation is suppressed 
significantly. On the other hand, comparative higher 
compression ratio, CR=21 leads to more NOx emission than that 
for CR=18. 



Crank Angle 



Fig. 6 Emission of NO x Vs Crank angle for different 
equivalence ratio and compression ratio 

4. Conclusion 

In the present study, a diesel engine has been investigated 
through numerical modelling focusing primarily on the 
influences of equivalence ratio and compression ratio on 
combustion and emission characteristics. In order to reduce the 
computational cost and time, a sector mesh of 45° angle with 
periodic boundary conditions applied at the periodic faces of the 
sector, is considered instead of using the whole engine geometry. 
CFD simulation are performed for a range of equivalence ratio 


varying from 0.6 to 1.2 and for three compression ratios namely, 
15:1, 18:1 and 21:1. For combustion performance, in-cylinder 
pressure, temperature and apparent heat release rate (AHRR) are 
analysed whereas emission of NOx, CO and unbumt 
hydrocarbon are also quantified to highlight the emission 
performance. Based on the modeling results, the major 
conclusions that can be drawn are followings: 

> Rapid combustion occurs for higher compression ratios 
and improvement in combustion characteristics with 
compression ratio for both lean and rich mixtures is 
observed. 

> Peak in-cylinder pressure and peak heat release nearer 
to TDC are achieved for compression ratio of 18:1 that 
could results in more engine torque. 

> For compression ratio beyond 16:1, rich mixture shows 
lower ignition delay compared to lean mixture. At 
lower compression ratio, in-cylinder temperature is not 
sufficiently high for atomization, vaporization, mixing 
of fuel with air, and preflame reactions to occur 
immediately after the fuel injection. 

> NOx emission in diesel engine increases due to higher 
pressure and temperature inside the cylinder associated 
with relatively higher compression ratio. Rich mixture 
results in more CO and unburnt hydrocarbon emission 
compared to lean mixture as result of incomplete 
combustion. 

> For fuel rich engine operation, a combustible mixture is 
readily formed in the combustion chamber and 
autoignition of fuel-air mixture occurs earlier. 
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